Package 'ss3sim'

Description

Add missing columns to each data frame in the list allowing for the use rbind() to create a single data frame. The code is based on rbind.fill from the plyr package.

Usage

add_colnames(dfs, bind = FALSE, fillwith = NA)

Arguments

Value

Depending on the input to bind you can either return the same structure, i.e., a list of data frames, or a data frame with all rows from each original data frame. Missing values will be filled with the entry in fillwith.

Author(s)

Kelli F. Johnson

Examples

x <- data.frame("a" = 1:10, "b" = 21:30)
y <- data.frame("a" = 11:15, "y" = letters[1:5])
alist <- ss3sim:::add_colnames(list(x, y), bind = FALSE)
adataframe <- ss3sim:::add_colnames(list(x, y), bind = TRUE)
# clean up
rm(x, y, alist, adataframe)

Description

Add short time varying parameter lines. At time of writing, this method will work for MG, selectivity, and catchability time varying, but not for SR

Usage

add_tv_parlines(string, tab, ctl_string, ss3.ctl)

Arguments

Value

A modified version of ss3.ctl (a vector of strings), containing the new parameter line

Description

Bias adjustment is performed to ensure that the only the most informative data available are used when estimating recruitment deviations. This process involves running the estimation method with the hessian both before this function and then running the estimation method again with the new values. Estimation files in the original folder will be deleted to ensure that convergence afterwards is not based on the input run. These files are permanently archived in the ⁠bias_[0-9]{2}⁠ folder within the directory. Only the three middle steps listed below are performed using calculate_bias, and the rest of the steps must be performed externally.

Usage

calculate_bias(dir, ctl_file_in)

Arguments

Details

• Estimate recruitment and the standard error about those estimates.

• Correct the estimates given their estimated uncertainty using a ramp.

• Save a new control file.

• Move original estimation files.

• Estimate model parameters.

Value

A list of bias adjustment parameters.

Author(s)

Kelli F. Johnson

Description

Given the sampling arguments that are specified in ..._params, e.g., index_params, calculate the super set of fleets, years, and data types that will be needed in the data file of expected values that is generated by the OM.

Usage

calculate_data_units(
  index_params = NULL,
  lcomp_params = NULL,
  agecomp_params = NULL,
  calcomp_params = NULL,
  mlacomp_params = NULL,
  wtatage_params = NULL
)

Arguments

Value

A list with the following three elements:

• fleets,

• years, and

• types.

Note

A superset by nature is larger than the individual sets used to create it (unless all sampling arguments are identical), so that the returned list will created some unnecessary combinations. This was done intentionally for simplicity but may be changed later.

Author(s)

Cole Monnahan

See Also

See further examples in clean_data and change_data

Examples

## Only one fleet
calculate_data_units(lcomp_params = list(fleets = 1, years = c(3, 4, 6)))
## Add new fleet
morefleets <- calculate_data_units(
  lcomp_params = list(fleets = 1, years = c(3, 4, 6)),
  agecomp_params = list(fleets = 2, years = 5)
)

## Add length or age if missing and conditional-age-at-length is included
test <- mapply(calculate_data_units,
  SIMPLIFY = FALSE,
  lcomp_params = list(NULL, list(fleets = 1, years = 1:10)),
  agecomp_params = list(NULL, NULL),
  MoreArgs = list(calcomp_params = list(fleets = 1, years = 1:10))
)

rm(test)

Description

Calculate the relative error (RE; $[EM - OM]/OM$) of parameters and derived quantities stored in a scalar or time series data frame generated by get_results_all().

Usage

calculate_re(dat, add = TRUE, EM = "em")

Arguments

Value

The default is to return a data frame structured the same as the input data frame, i.e., dat, but with additional columns, where '_re' is appended to the base string of the column name. All NAN and Inf values are returned as NA values, typically because you cannot divide by zero. Irrelevant columns, i.e., columns of entirely zero of NA are removed prior to returning the data frame.

Author(s)

Sean Anderson and Cole Monnahan

See Also

get_results_all(), get_results_scenario()

Examples

# Example with built in package data:
data("ts_dat", package = "ss3sim")
data("scalar_dat", package = "ss3sim")
head(calculate_re(ts_dat))
head(calculate_re(ts_dat, add = FALSE))
head(calculate_re(scalar_dat, add = FALSE))
rm("ts_dat", "scalar_dat")

Description

Change catch in the data so at least all combinations of fleet, season, and year, needed for catch are available. Equilibrium years are generated if there are equilibrium parameters in the control list.

Usage

change_catch(dat_list, ctl_list)

Arguments

Value

A modified Stock Synthesis data file as a list in R.

Author(s)

Kathryn L. Doering

See Also

change_f() changes the fishing mortality, $F$, parameters using the control file, but these $F$ values will only be implemented for years with corresponding entries in the Stock Synthesis data file. Thus, this function must be implemented after change_f().

Description

Change the composition data in a Stock Synthesis data list object or file to include rows of data that are desired. Typically, this will be an operating model (OM) because only dummy-data observations are used here, i.e., all compositions are set to a value of one. Creating these dummy observations is helpful before running your OM because it will facilitate the creation of observed values for each desired combination.

Usage

change_comp(dat_list, type = c("len", "age", "cal"), paramlist, nsex = 1, bins)

Arguments

Examples

# todo: remove this example when testing is complete
## Not run: 
change_comp(
  dat_list = dat, type = "len",
  paramlist = scenariol[[1]][[c("lcomp_params", "agecomp_params")]]
)

## End(Not run)

Description

Change the bins of a data frame object from a dat_list, filling in the columns with ones.

Usage

change_dat_bin(object, bins)

Arguments

Value

A modified data frame where columns holding old composition data are removed in their entirety and new columns of ones are filled for each value in bins.

Description

Alter the structure of data that is available from a Stock Synthesis operating model (OM), which in turn leads to changes in the output and ability to sample data after running the model.

Usage

change_data(
  dat_list,
  outfile = NULL,
  fleets,
  years,
  types = c("len", "age", "cal", "mla", "mwa"),
  age_bins = NULL,
  len_bins = NULL,
  pop_binwidth = NULL,
  pop_minimum_size = NULL,
  pop_maximum_size = NULL,
  lcomp_constant = NULL,
  tail_compression = NULL,
  nsex = 1
)

Arguments

Details

change_data() is called internally within ss3sim, but it can be used to manipulate data or to prepare a new OM for use in a simulation. Original data is removed and dummy data is added to the Stock Synthesis .dat object. The dummy data expands the data structure to provide information for all years and fleets, potentially adding many rows of data.

Currently, .dat files with multiple sexes cannot be manipulated with change_data().

The robustification constant is added to both the observed and expected proportions of length composition data, before being normalized internally. It is designed to help stabilize the model, but is unclear how and when to use it for optimal effect. The same value is used for all length data.

Value

An invisible data list, and a file is written to the disk if an entry other than the default of NULL is provided for outfile.

Author(s)

Cole Monnahan, Ian G. Taylor, Sean Anderson, Kelli F. Johnson

See Also

See clean_data() for a counter function.

Other change functions: change_e(), change_em_binning(), change_f(), change_o(), change_retro(), change_tv()

Description

Takes Stock Synthesis .ctl and forecast.ss files, along with a list structure which houses the data file as read in by r4ss::SS_readdat() and changes which parameters are estimated, how natural mortality is estimated, and if forecasts are performed. The function can be called by itself or within run_ss3sim() to alter an estimation model .ctl file.

Usage

change_e(
  ctl_file_in = "em.ctl",
  ctl_file_out = "em.ctl",
  dat_list = NULL,
  for_file_in = "forecasts.ss",
  par_name = NULL,
  par_int = "NA",
  par_phase = "NA",
  forecast_num = 0,
  verbose = FALSE
)

Arguments

Value

Altered versions of Stock Synthesis .ctl and forecast.ss files are written to the disk and the altered dat_list is returned invisibly.

Author(s)

Kelli F. Johnson

See Also

Other change functions: change_data(), change_em_binning(), change_f(), change_o(), change_retro(), change_tv()

Examples

d <- system.file("extdata", "models", "cod-om", package = "ss3sim")
change_e(
  ctl_file_in = file.path(d, "codOM.ctl"),
  ctl_file_out = file.path(tempdir(), "change_e.ctl"),
  dat_list = codomdat,
  for_file_in = file.path(d, "forecast.ss"),
  par_name = c("_steep", "Size_DblN_peak_Fishery(1)"),
  par_int = c(0.3, 40), par_phase = c(3, 2),
  forecast_num = 0
)
# clean up the temporary files
file.remove(file.path(tempdir(), "change_e.ctl"))

Description

change_em_binning() alters the bin structure for the population and length-composition data in a Stock Synthesis estimation model (EM). The original length-composition data from the EM .dat is changed according to the user's specification. If the data file also contains conditional age-at-length data, then these data will be re-binned as well.

Usage

change_em_binning(
  dat_list,
  outfile = NULL,
  bin_vector,
  lbin_method = NULL,
  pop_binwidth = NULL,
  pop_minimum_size = NULL,
  pop_maximum_size = NULL
)

Arguments

Author(s)

Kotaro Ono (length-composition rebinning), Sean Anderson (conditional age-at-length rebinning)

See Also

Other change functions: change_data(), change_e(), change_f(), change_o(), change_retro(), change_tv()

Examples

# Note that typically this function is used with estimation models in ss3sim,
# but it is used with an operating model data file in the following examples.
f <- system.file("extdata", "models", "cod-om", "codOM.dat", package = "ss3sim")
d <- r4ss::SS_readdat(f, verbose = FALSE)

# An example with lbin_method = 1
l1 <- change_em_binning(d,
  outfile = NULL, lbin_method = 1,
  bin_vector = seq(20, 152, by = 4)
)
l1$lbin_vector
head(l1$lencomp)

# An example with lbin_method = 2
new_bin_vec <- seq(min(d$lbin_vector), max(d$lbin_vector), by = 4)
# add the max value if necessary.
if (new_bin_vec[length(new_bin_vec)] != d$lbin_vector[length(d$lbin_vector)]) {
  new_bin_vec <- c(
    new_bin_vec,
    d$lbin_vector[length(d$lbin_vector)]
  )
}
pop_bin_input <- 5
pop_min_size_input <- min(d$lbin_vector_pop) - 1
pop_max_size_input <- max(d$lbin_vector_pop) + 5
lbin_vec_pop <- seq(pop_min_size_input,
  pop_max_size_input,
  length.out = (pop_max_size_input - pop_min_size_input) /
    pop_bin_input + 1
)
l2 <- change_em_binning(
  dat_list = d,
  bin_vector = new_bin_vec,
  lbin_method = 2,
  # Note: need more inputs with lbin_method = 2
  pop_binwidth = pop_bin_input,
  pop_minimum_size = pop_min_size_input,
  pop_maximum_size = pop_max_size_input
)
l2$lbin_method
# note bin width is now the same as the input
pop_bin_input
l2$binwidth
# note the minimum size has changed based on the input:
pop_min_size_input
l2$minimum_size
# so has max
l2$maximum_size
l2$lbin_vector
# other modified components:
l2$lbin_vector_pop
head(l2$lencomp)

Description

Replace or input a time series of fishing mortality, $F$, values into a Stock Synthesis control file. In Stock Synthesis, inserting $F$ values in this manner, relies on the assumption that $F$ operates continuously throughout the year and and the process operates jointly with natural mortality (Baranov 1918; Branch 2009). The documentation for Stock Synthesis also describes this process as $F$ method == 2, where $F$ is continuous and modeled using full parameters.

Usage

change_f(years, fleets, fvals, seasons = 1, ses = 0.005, ctl_list)

Arguments

Details

The argument years is the only argument that must be a vector or a list of vectors. Other arguments can be specified using a single scalar value that will be repeated for all fisheries in all years. If the input argument needs to be different for any year or fishery, the argument must be a list with vectors for each fishery, where each vector is the same length as the vectors within the years argument. Although, both the years and other input arguments can be specified using a single vector if the length of fleets is just one or a vector of values is specified for fleets because all of these vectors will just be combined into a single data frame. Where it gets complicated is when there are multiple fleet and year combinations, then it is best to just use the list structure common to other functions within ss3sim.

change_f() overrides any $F$ values that are in the supplied control file with the newly specified values, i.e., fvals. Users do not need to specify values for years in which there will be zero fishing because Stock Synthesis will automatically set them to zero when running the operating model. Using the control file rather than the par file to manipulate the operating model requires a few other files within the operating model folder to be set up in a particular manner. That is, (1) the starter file must be set up to read parameters from the control file rather than the par file and (2) the data file must have a dummy catch entry for every year, fishery combination that will be specified in the control file. If a year, fishery combination is specified in the control file and not present in the data file, then the entry in the control file will be ignored. ss3sim automatically corrects for this using ss3sim_base() by specifying a row for every year and fleet using change_catch().

Value

Modified Stock Synthesis control file list.

Author(s)

Kelli F. Johnson

See Also

See r4ss::SS_readctl() and r4ss::SS_writectl() for how to supply ctl_list and how to write the file back to the disk once you are done manipulating the list object.

Other change functions: change_data(), change_e(), change_em_binning(), change_o(), change_retro(), change_tv()

Examples

dat <- r4ss::SS_readdat(
  system.file("extdata", "models", "cod-om", "codOM.dat", package = "ss3sim"),
  verbose = FALSE
)
ctl <- r4ss::SS_readctl(
  system.file("extdata", "models", "cod-om", "codOM.ctl", package = "ss3sim"),
  verbose = FALSE, use_datlist = TRUE, datlist = dat
)
# Using original vector-style inputs
newctl <- change_f(years = 1:50, fleets = 1, fvals = 0.2, ctl_list = ctl)
# Using list-style inputs for when there are multiple fisheries
newctl <- change_f(
  years = list(1:5, 1:10), fleets = 3:4,
  fvals = list(rep(0.1, 5), rep(0.2, 10)), ctl_list = ctl
)
rm(dat, ctl, newctl)

Description

This function replaces the robustification value for length-composition data in a .dat file that was read in using r4ss::SS_readdat() with those specified in lcomp_constant. It then writes a new file with name outfile into the working directory.

Usage

change_lcomp_constant(lcomp_constant, dat_list, outfile = NULL)

Arguments

Details

The robustification constant is added to both the observed and expected proportions of length composition data, before being normalized internally. It is designed to help stabilize the model, but is unclear how and when to use it for optimal effect. The same value is used for all length data.

Value

A modified Stock Synthesis .dat file, and that file returned invisibly (for testing) as a vector of character lines.

Author(s)

Cole Monnahan

Description

change_o takes a Stock Synthesis .ctl file and implements parameter value changes that are NOT time varying. change_o is specifically set up to work with an operating model .ctl file.

Usage

change_o(
  change_o_list,
  ctl_file_in = "control.ss_new",
  ctl_file_out = "om.ctl",
  par_name = NULL,
  par_int = NULL,
  verbose = FALSE
)

Arguments

Value

The function creates modified versions of the .ctl files. The function also returns change_o_list invisibly.

Specifying the change_o_list

Parameters initial values will change according to the values passed to change_o_list. Each parameter should have a single value specified.

Parameter names must be unique and match the full parameter name in the .ctl file.

Passing arguments to change_o() through run_ss3sim()

(1) add a column called co.par_name to the simdf that specifies which parameters you want to change in the OM, each element of this vector needs to be wrapped in quotations to be later evaluated, e.g., 'c("SR_BH_steep","SR_sigmaR")' represents a single entry; and (2) add an additional column called co.par_int to the simdf that specifies INIT values for each parameter in the previous column, e.g., "c(0.6, 1.0)", if there is more than one value, the vector needs to be wrapped in quotations to be evaluated later.

Author(s)

Kathryn L. Doering

See Also

Other change functions: change_data(), change_e(), change_em_binning(), change_f(), change_retro(), change_tv()

Description

The population length bins in Stock Synthesis structure size data and empirical weight-at-age data. change_pop_bin changes the data file to contain specifications to create a vector (length-bin method of 2) rather than the actual bins from the length data (length-bin method of 1) or an actual vector (length-bin method of 3).

Usage

change_pop_bin(
  dat_list,
  binwidth = NULL,
  minimum_size = NULL,
  maximum_size = NULL,
  maximum_age = NULL
)

Arguments

Details

The only required argument is dat_list and the remaining arguments default to a value of NULL, which leads to the data file not being changed.

Value

A modified Stock Synthesis data file in list form. The list is only returned if it is assigned to an object.

Description

Manipulates the control list to simultaneously add and remove elements related to time series data on trends.

Usage

change_q(
  string_add = NULL,
  string_remove = NULL,
  ctl_list,
  dat_list = lifecycle::deprecated(),
  ctl_file_in = lifecycle::deprecated(),
  dat_file_in = lifecycle::deprecated(),
  ctl_file_out = lifecycle::deprecated(),
  overwrite = lifecycle::deprecated(),
  verbose = lifecycle::deprecated()
)

Arguments

Details

Catchability, $q$, represents the proportionality constant between data on trends and estimated population abundance. Thus a survey thought to encapsulate the entire population, e.g., an acoustic survey of the entire area, will have $q = 1$. In Stock Synthesis, environmental time series are modelled similarly to a survey or catch-per-unit-effort time series and thus will also have a catchability term. Readers interested in the complete range of functionality should see the catchability section of the Stock Synthesis user manual. change_q() has limited functionality relative to what is available in Stock Synthesis. For example, change_q() cannot add parameters for additional variance. Though it will remove additional variance parameters for fleets that no longer have survey data. Additionally, the float term is not used within ss3sim and is set to zero.

Value

A modified Stock Synthesis control list.

Author(s)

Kelli F. Johnson

See Also

• check_q() determines which fleets should removed or added.

• r4ss::SS_readctl() reads in the control file passed to ctl_list.

• find_position() allows ⁠string_*⁠ to use strings or integers.

Examples

removedfleet1 <- change_q(string_remove = 1, ctl_list = codomctl)
removedfleet2 <- change_q(string_remove = 2, ctl_list = codomctl)
removedfleets <- change_q(
  string_remove = c("Fishery", 2),
  ctl_list = codomctl
)
testthat::expect_null(removedfleets[["Q_options"]])
newctl <- codomctl
newctl[["fleetnames"]] <- c(newctl[["fleetnames"]], "testfleet")
newctl[["Nfleets"]] <- length(newctl[["fleetnames"]])
newctl <- change_q(string_add = "testfleet", ctl_list = newctl)
testthat::expect_equal(newctl[["Q_options"]][, "fleet"], 1:3)

Description

This function replaces the recruitment deviations in the control file of a Stock Synthesis model with those specified in the argument recdevs. The new control file is then written to the disk if ctl_file_out is specified. It is imperative that the path provided in ctl_file_in be to a ss_new file so change_rec_devs can properly determine where to place the recruitment deviations in the control file.

Usage

change_rec_devs(recdevs, ctl_file_in, ctl_file_out = "control_recruitment.ss")

Arguments

Value

A modified Stock Synthesis control file.

Author(s)

Kelli F. Johnson

Examples

d <- system.file(file.path("extdata", "models"), package = "ss3sim")
change_rec_devs(
  recdevs = rlnorm(101),
  ctl_file_in = file.path(d, "cod-om", "codOM.ctl"),
  ctl_file_out = file.path(tempdir(), "control_recdevs.ss")
)
# Change the recruitment deviations in years 2:11
change_rec_devs(
  recdevs = stats::setNames(rlnorm(10), 2:11),
  ctl_file_in = file.path(d, "cod-om", "codOM.ctl"),
  ctl_file_out = file.path(tempdir(), "control_recdevsInitial.ss")
)
lapply(
  X = dir(tempdir(), pattern = "control_.+ss", full.names = TRUE),
  FUN = unlink
)

Description

Change start year main recruitment deviations in control file

Usage

change_recyear(ctl_list, main)

Arguments

Value

A r4ss::SS_readctl() list with an augmented start year of the recruitment deviations in the main period.

Author(s)

Kelli F. Johnson

Description

A retrospective analysis tests the effect of peeling back the number of operating model years observable to the estimation model. This function alters the Stock Synthesis starter file to run a retrospective analysis.

Usage

change_retro(
  str_file_in = "starter.ss",
  str_file_out = "starter.ss",
  retro_yr = 0
)

Arguments

Details

Note that the starter file is set up to run a single retrospective run. Therefore, if you would like to run retrospective analyses for, say, 0, 1, 2, 3, 4, and 5 years, you will need to use this function to adjust the starter file 6 separate times.

Value

A modified Stock Synthesis starter file.

Author(s)

Sean C. Anderson

See Also

Other change functions: change_data(), change_e(), change_em_binning(), change_f(), change_o(), change_tv()

Examples

# Create a temporary folder for the output:
temp_path <- file.path(tempdir(), "ss3sim-retro-example")
dir.create(temp_path, showWarnings = FALSE)

# Locate the package data:
starterfile <- system.file("extdata", "models", "cod-om",
  "starter.ss",
  package = "ss3sim"
)

# No retrospective analysis:
change_retro(starterfile, paste0(temp_path, "/retro-0-starter.ss"),
  retro_yr = 0
)

# A retrospective analysis of 5 years:
change_retro(starterfile, paste0(temp_path, "/retro-5-starter.ss"),
  retro_yr = -5
)

Description

Change start year of the data file

Usage

change_startyear(dat_list, firstyear = NULL)

Arguments

Value

A r4ss::SS_readdat() list with an augmented start year.

Author(s)

Kelli F. Johnson

Description

This function replaces the tail compression value for length-composition data in a .dat file that was read in using r4ss::SS_readdat() with those specified in tail_compression. It then writes a new file with name dat_file_out into the working directory.

Usage

change_tail_compression(tail_compression, dat_list, outfile = NULL)

Arguments

Value

A modified Stock Synthesis .dat file is returned invisibly.

Author(s)

Cole Monnahan

Description

change_tv takes Stock Synthesis .ctl, .par, and .dat files and implements time-varying parameters using environmental variables. change_tv is specifically set up to work with an operating model .ctl file.

Usage

change_tv(
  change_tv_list,
  ctl_file_in = "control.ss_new",
  ctl_file_out = "om.ctl",
  dat_file_in = "ss3.dat",
  dat_file_out = "ss3.dat"
)

Arguments

Details

Although there are three ways to implement time-varying parameters within Stock Synthesis, ss3sim and change_tv only use the environmental variable option. Within Stock Synthesis, time-varying parameters work on an annual time-step. Thus, for models with multiple seasons, the time-varying parameters will remain constant for the entire year.

The ctl_file_in argument needs to be a .ss_new file because the documentation in .ss_new files are automated and standardized. This function takes advantage of the standard documentation the .ss_new files to determine which lines to manipulate and where to add code in the .ctl, .par, and .dat files, code that is necessary to implement time-varying parameters.

ss3sim uses annual recruitment deviations and may not work with a model that ties recruitment deviations to environmental covariates. If you need to compare the environment to annual recruitment deviations, the preferred option is to transform the environmental variable into an age 0 pre-recruit survey. See page 55 of the Stock Synthesis version 3.24f manual for more information.

Value

The function creates modified versions of the .ctl and .dat files if ctl_file_out and dat_file_out are not NULL. The function also returns a list of the modified .ctl and .dat R objects invisibly.

Specifying the change_tv_list

Parameters will change to vary with time according to the vectors of deviations passed to change_tv_list. Vectors of deviations, also referred to as environmental data, must have a length equal to endyr-startyr+1, where endyr and startyr are specified the .dat file. Specify years without deviations as zero.

Parameter names must be unique and match the full parameter name in the .ctl file. Names for stock recruit parameters must contain "devs", "R0", or "steep", and only one stock recruit parameter can be time-varying per model.

This feature will include an additive functional linkage between environmental data and the parameter where the link parameter is fixed at a value of one and the par value is specified in the .par file: $par'[y] = par + link * env[y]$.

For catchability ($q$) the additive functional linkage is implemented on the log scale: $ln(q'[y]) = ln(q) + link * env[y]$

Author(s)

Kotaro Ono, Carey McGilliard, Kelli F. Johnson, and Kathryn L. Doering

See Also

Other change functions: change_data(), change_e(), change_em_binning(), change_f(), change_o(), change_retro()

Examples

## Not run: 
# Create a temporary folder for the output and set the working directory:
temp_path <- file.path(tempdir(), "ss3sim-tv-example")
dir.create(temp_path, showWarnings = FALSE)
wd <- getwd()
setwd(temp_path)
on.exit(setwd(wd), add = TRUE)

d <- system.file("extdata", package = "ss3sim")
om <- file.path(d, "models", "cod-om")
dir.create("cod-om")
file.copy(om, ".", recursive = TRUE)
setwd("cod-om")

change_tv(
  change_tv_list =
    list(
      "NatM_uniform_Fem_GP_1" = c(rep(0, 20), rep(.1, 80)),
      "SR_BH_steep" = stats::rnorm(100, 0, 0.05)
    ),
  ctl_file_in = "codOM.ctl",
  ctl_file_out = "example.ctl",
  dat_file_in = "codOM.dat",
  dat_file_out = "example.dat"
)

# Clean up:
unlink("cod-om", recursive = TRUE)

## End(Not run)

Description

Keep all of the data in the model but change the years that are estimated in the model. First year of the model will be first year of non-zero catch. Main recruitment period starts 1/2 generation time before first year of compositional data included in the model. Late recruitment is the last year of the model by default and cannot be modified using this function, neither can early recruitment, which starts in year 1.

Usage

change_year(dat_list, ctl_list)

Arguments

Description

Check that the Stock Synthesis data file looks correct

Usage

check_data(x)

Arguments

Description

Check that the param list inputs have correct structure and range given an associated data file.

Usage

check_data_str_range(all_params, dat_list)

Arguments

Description

Calculate the length of all input arguments to see if they are equal. Entries that are NULL, and thus, have a length of zero are ignored. An optional trigger to stop() is provided with a tailored error message.

Usage

check_eqlength(..., keepgoing = FALSE)

Arguments

Value

TRUE or FALSE depending on the result of the test. Nothing is returned if the stop function is invoked.

Author(s)

Kelli F. Johnson

Examples

Description

Ensure that the forecast.ss file is configured for use in ss3sim.

Usage

check_forecast(for_list)

Arguments

Details

• fish at $F_{MSY}$

• use relative benchmark years (i.e., Bmark_years)

• use relative years for fishing specifications, i.e., Fcast_years

Value

A an augmented list object, as returned by r4ss::SS_readforecast(), is invisibly returned.

Author(s)

Kelli F. Johnson

Description

Check a Stock Synthesis control file to determine if the desired fleets have q parameters set up.

Usage

check_q(ctl_list, Nfleets = lifecycle::deprecated(), desiredfleets)

Arguments

Value

A list with two vectors, add and remove, specifying which fleets to add and which to remove from the control file.

See Also

change_q() for actually adding or removing the fleets.

Examples

# Keep just the fishery
stopifnot(check_q(ctl_list = codomctl, desiredfleets = 1)[["remove"]] == 2)
# All elements of the returned list should be NULL
# because the model only has two \eqn{q} parameters
stopifnot(all(mapply(is.null, check_q(codomctl, desiredfleets = 1:2))))
# Fleet 3 is not present
stopifnot(check_q(codomctl, desiredfleets = 1:3)[["add"]] == 3)
stopifnot(check_q(codomctl, desiredfleets = 2:3)[["remove"]] == 1)

Description

This prepares a .dat file to be used by an estimation method, whereas before it may have had leftover data from sampling purposes.

Usage

clean_data(
  dat_list,
  lcomp_params = NULL,
  agecomp_params = NULL,
  calcomp_params = NULL,
  mlacomp_params = NULL,
  verbose = FALSE
)

Arguments

Value

An invisible cleaned data list as an object.

Note

This function does not write the result to file.

Author(s)

Cole Monnahan

See Also

calculate_data_units()

Other sampling functions: sample_agecomp(), sample_calcomp(), sample_catch(), sample_discard(), sample_index(), sample_lcomp(), sample_mlacomp(), sample_wtatage()

Description

A list of controls returned from r4ss::SS_readctl() for the North Sea cod operating model. The input file is stored in extdata/models.

Usage

codemctl

Format

A list with many items, some of which are highlighted below:

fleetnames

a vector of names for the fleets

MG_parms

natural mortality and growth parameters

SR_parms

stock-recruitment relationship parameters

...

Source

North Sea cod (Gadus morhua; Richard D. Methot, Jr., NMFS, NOAA, pers. comm.)

See Also

• r4ss::SS_readctl()

• r4ss::SS_readdat()

Examples

data("codemctl", package = "ss3sim")

Description

A list of controls returned from r4ss::SS_readctl() for the North Sea cod operating model. The input file is stored in extdata/models.

Usage

codomctl

Format

A list with many items, some of which are highlighted below:

fleetnames

a vector of names for the fleets

MG_parms

natural mortality and growth parameters

SR_parms

stock-recruitment relationship parameters

...

Source

North Sea cod (Gadus morhua; Richard D. Methot, Jr., NMFS, NOAA, pers. comm.)

See Also

• r4ss::SS_readctl()

• r4ss::SS_readdat()

Examples

data("codomctl", package = "ss3sim")

Description

A list of data returned from r4ss::SS_readdat() for the North Sea cod operating model. The input file is stored in extdata/models.

Usage

codomdat

Format

A list with many items, some of which are highlighted below:

catch

data frame of catches by year, fleet, and season

CPUE

catch-per-unit-effort data

lencomp

length-composition data

...

Source

North Sea cod (Gadus morhua; Richard D. Methot, Jr., NMFS, NOAA, pers. comm.)

See Also

r4ss::SS_readdat()

Examples

data("codomdat", package = "ss3sim")

Description

This function exists for back compatibility. Note that this will only work if the column model_run has only the strings"om" or "em".

Usage

convert_to_wide(lng)

Arguments

Value

A wide dataframe (separate columns for em and om results)

Author(s)

Kathryn L. Doering

Examples

## Not run: 
scalar <- utils::read.csv("ss3sim_scalar.csv")
scalar_wide <- convert_to_wide(scalar)

ts <- utils::read.csv("ss3sim_ts.csv")
ts_wide <- convert_to_wide(scalar)

## End(Not run)

Description

Copy the OM or EM into a scenario directory

Usage

copy_ss3models(model_dir, scenario, iteration = 1, type = c("om", "em"))

Arguments

Value

An invisible boolean for whether that iteration already existed.

Author(s)

Sean C. Anderson, Kelli F. Johnson

Examples

# Locate the package data:
om_folder <- system.file(
  "extdata", "models", "cod-om",
  package = "ss3sim"
)

# Copy the operating model:
copy_ss3models(
  model_dir = om_folder,
  scenario = "D0-F0-testing"
)
# Now look at your working directory in your file system

# Copy the EM
copy_ss3models(
  model_dir = om_folder,
  type = "em",
  scenario = "D1-F0-testing"
)
# Scenario argument affects the folder names.

# Clean up:
unlink("D0-F0-testing", recursive = TRUE)
unlink("D1-F0-testing", recursive = TRUE)

Description

Create estimation model (EM) files from operating model (OM) files. By making small changes to the OM rather than having two sets of files, less files need to be maintained. Differences between the OM and EM are mainly related to how the OM takes input fishing mortality values rather than absolute catches.

Usage

create_em(
  dir_in = system.file("extdata", "models", "cod-om", package = "ss3sim"),
  dir_out = file.path(getwd(), "new-em")
)

Arguments

Value

Nothing is returned, but three files are saved to the disk in the specified folder that may also be new.

Control file

Most changes to the EM control file relate to recruitment and fishing. The phase in which recruitment deviations are estimated is checked to ensure that it is positive. Though, this might be unnecessary because the OM file can have negative or positive phases. Thus, users are encouraged to just set the phase in which recruitment is estimated in the OM at the value that they would like to use in the EM. Additional changes are made to the bias adjustment procedure based on the biology of the stock.

The F_Method is set to 3 to allow the model to estimate fishing mortality based on catches in the data file. Users might want to adjust the maximum fishing mortality based on their scenarios.

Data file

No data file is needed for the EM. The data_expval.ss file produced when executing the OM contains the expected values of the OM population dynamics. ss3sim provides three functions which carry out the random sampling process and generate .dat files to be used in the EM. See the Introduction vignette vignette("introduction", package = "ss3sim") for more details.

Forecast file

Nothing is changed in the forecast file from the OM.

Starter file

The names of the data and control files are specified and the maximum phase for estimation is set to 100.

Author(s)

Kelli F. Johnson

Examples

create_em()
# The necessary files are in the following folder
dir(file.path(getwd(), "new-em"))
# Clean up your directory
unlink(file.path(getwd(), "new-em"), recursive = TRUE)

Description

Generate and save, if outfile is provided, the ss3sim logo using the built-in data.

Usage

create_logo(outfile = NULL)

Arguments

Value

A png file or a graphics device with the logo used for the ss3sim project.

Author(s)

Kelli F. Johnson

Examples

ss3sim:::create_logo()
grDevices::dev.off()

Description

Used internally by the plotting functions to create faceting formulas.

Usage

facet_form(horiz = NULL, horiz2 = NULL, vert = NULL, vert2 = NULL)

Arguments

Value

A formula which can be used in ggplot2::facet_grid() or NULL if all arguments are NULL.

Author(s)

Cole Monnahan

Description

Function that fills in matrix across rows of wtatage data by interpolation Missing Rows are then backfilled

Usage

fill_across(mat, minYear, maxYear)

Arguments

Author(s)

Peter Kuriyama and Allan Hicks

See Also

• sample_lcomp()

• sample_agecomp()

• sample_wtatage()

Description

Find the position of each desired value, i.e., x, in a vector of strings. Builds on match() by allowing x to be a combination of strings to be matched and known positions.

Usage

find_position(x, table)

Arguments

Value

An integer vector indicating the positions of x in table. Same as match(), if x[i] is found to be equal to table[j], then the value returned in the i-th position of the integer vector is j. The smallest value of j, i.e., the first match, is always returned.

find_position() differs from match() in three ways. First, values of x that are not found are removed. Thus, the length of the integer vector has the potential to be shorter than the length of x. Second, x can contain a mix of integer positions that are already known and strings to be found. Third, table cannot include integers that do not match their position. See the specifications for table for more details.

Author(s)

Kelli F. Johnson

See Also

See match() for a more formal version of find_position() that returns an integer vector the same length as x.

Examples

# Standard use
find_position(c("sad", 1), c("happy", "sad"))
# Incorrect use
find_position(c("sad", 2), c("happy", "sad", "2"))

Description

Get Stock Synthesis binary/executable location

Usage

get_bin(bin_name = "ss3")

Arguments

Value

A string providing the full path to a Stock Synthesis binary. If using the GitHub version of ss3sim, this will be an internal binary. Otherwise, get_bin() searches for a version of the binary in your path. See the ss3sim vignette fore more information.

Author(s)

Sean C. Anderson

Examples

## Not run: 
get_bin()

## End(Not run)

Description

Extract the summary of fits to composition data, where the sections are structured similarly for each type of data in the report file.

Usage

get_compfit(report.file, name)

Arguments

Description

This function returns the location of one of the built-in model configurations.

Usage

get_model_folder(folder_name)

Arguments

Value

A character object showing the location of the appropriate model configuration folder in the package extdata folder.

Examples

get_model_folder("cod-em")

Description

Names of the available NLL components will depend on the version of the model. Names are native to the estimation framework and all available components are extracted.

Usage

get_nll_components(report.file)

Arguments

Value

A vector of named numeric values, where "NLL_" is appended to the names in the report.file.

Author(s)

Merrill Rudd

Description

This function returns a set of pseudo-random recruitment deviations based on an iteration number. Given the same iteration number the function will return the same recruitment deviations. The deviations are standard normal. I.e., they have a mean of 0 and a standard deviation of 1.

Usage

get_recdevs(iteration, n, seed = 21)

Arguments

Value

A vector of standard normal recruitment deviations.

Examples

get_recdevs(1, 10)
get_recdevs(1, 10)
get_recdevs(2, 10)

Description

This high level function extracts results from Stock Synthesis model runs. Give it a directory which contains directories for different "scenario" runs, within which are iterations. It writes two data.frames to file: one for single scalar values (e.g., MSY) and a second that contains output for each year of the same model (timeseries, e.g., biomass(year)). These can always be joined later.

Usage

get_results_all(
  directory = getwd(),
  overwrite_files = FALSE,
  user_scenarios = NULL,
  type = c("long", "wide"),
  filename_prefix = "ss3sim"
)

Arguments

Value

Returns a list of 3 dataframes: scalar, ts, and dq. Creates two .csv files in the current working directory, where the names of those files are based on filename_prefix and the default leads to the following: ss3sim_ts.csv and ss3sim_scalar.csv.

Author(s)

Cole Monnahan, Merrill Rudd, Kathryn L. Doering

See Also

Other get-results: get_results_derived(), get_results_scalar(), get_results_scenario(), get_results_timeseries()

Description

Extract time series from an r4ss::SS_output() list from a model run. Returns a data.frame of the results for spawning stock biomass (SSB), recruitment, forecasts, and effort by year.

Usage

get_results_derived(report.file)

Arguments

Author(s)

Kelli F. Johnson

See Also

Other get-results: get_results_all(), get_results_scalar(), get_results_scenario(), get_results_timeseries()

Description

Get results for 1 iteration

Usage

get_results_iter(dir_1_iter = NULL, mod_dirs = NULL, iter_name = NULL)

Arguments

Value

A list of 3 data frames called scalar, timeseries, and derived (for derived quantities). These lists contain information for multiple model runs (estimation models and operating models) for 1 iteration.

Author(s)

Kathryn L. Doering

Description

Get results for 1 model run

Usage

get_results_mod(dir = getwd(), is_EM = NULL, is_OM = NULL)

Arguments

Value

A list of 3 data frames called scalar, timeseries, and derived (for derived quantities). These data frames contain results for 1 model run.

Author(s)

Kathryn L. Doering

Description

Extract scalar quantities from an r4ss::SS_output() list from a model run. Returns a data.frame of the results (a single row) which can be rbinded later.

Usage

get_results_scalar(report.file)

Arguments

Author(s)

Cole Monnahan; Merrill Rudd

See Also

Other get-results: get_results_all(), get_results_derived(), get_results_scenario(), get_results_timeseries()

Description

Extract results from all iterations inside a supplied scenario folder. The function writes the following .csv files to the scenario folder

1. scalar metrics with one value per iteration (e.g. $R_0$, $h$),

2. a timeseries data ('ts') which contains multiple values per iteration (e.g. $SSB_y$ for a range of years $y$), and

3. residuals on the log scale from the surveys across all iterations; this functionality is currently disabled and not tested.

Usage

get_results_scenario(scenario, directory = getwd(), overwrite_files = FALSE)

Arguments

Author(s)

Cole Monnahan and Kathryn L. Doering

See Also

get_results_all() loops through these .csv files and combines them together into a single "final" dataframe.

Other get-results: get_results_all(), get_results_derived(), get_results_scalar(), get_results_timeseries()

Description

Extract and return time series from an r4ss::SS_output() list, that is read in from the estimation method of a single iteration. The main time series information is included but no information about the uncertainty of those measurements is available. See the derived quantities for uncertainty.

Usage

get_results_timeseries(report.file)

Arguments

Details

Information about both season and area are included in the data frame. For values that have no associated season or area, i.e., are summary values over all areas and seasons, the values are repeated for each area/season combination within a given year. For example, the recruitment deviation is for all areas and is thus repeated in each row across areas for a given year.

Value

A data frame with the following columns:

• year

• Area

• Seas

• Bio_smry

• SpawnBio

• Recruit_0

• retainB_0-9+

• retainN_0-9+

• deadB_0-9+

• deadN_0-9+

• F_0-9+

• SPRratio

• rec_dev

• raw_rec_dev

Author(s)

Cole Monnahan

See Also

Other get-results: get_results_all(), get_results_derived(), get_results_scalar(), get_results_scenario()

Description

Find scenario directories, where it is known if a directory was derived from ss3sim and contains iterations of the operating and estimating models if there are directories with numeric names that contain om and/or em directories, i.e., "1", "2", "3", ..., "100".

Usage

get_scenarios(directory = getwd(), full = FALSE)

Arguments

Value

A character vector of names of directories that contain output from ss3sim. Full paths are only provided if full = TRUE.

Author(s)

Merrill Rudd

Description

Use the name of the operating model to open the ctl file and obtain the INIT value for sigmaR (recruitment deviations sigma)

Usage

get_sigmar(om)

Arguments

Author(s)

Kelli F. Johnson

Description

Use the presence of files generated by Stock Synthesis to determine if the run was successful or not. There are two levels of success that must be determined if the run was meant to include estimating the hessian.

Usage

get_success(dir)

Arguments

Value

A named vector with two values, ran and hess, (1) a zero or one for the presence of the Report.sso file and (2) a zero or one for the presence of a positive definite hessian matrix.

Author(s)

Kelli F. Johnson

Description

Bind together list of list components with the same name

Usage

make_df(list_name, list_df)

Arguments

Value

A dataframe

Author(s)

Kathryn L. Doering

Description

Generate boxplots using ggplot2::ggplot() to visualize outliers and central tendencies.

Usage

plot_boxplot(
  data,
  x,
  y,
  horiz = NULL,
  horiz2 = NULL,
  vert = NULL,
  vert2 = NULL,
  relative.error = FALSE,
  axes.free = TRUE,
  print = TRUE,
  fill = NA
)

Arguments

Details

Median, hinges, and whiskers as well as outliers are displayed to summarize the data. The lower and upper hinges are the first and third quantiles (i.e., 25th and 75th percentiles). The upper and lower whiskers are 1.5*inner-quartile range, i.e., the distance between the first and third quartiles. Outliers are those points that lie beyond the whiskers. These explanations are detailed in ggplot2::geom_boxplot().

Values of NA are removed prior to plotting such that the typical error message from ggplot2::ggplot() is not printed to the screen.

The ss3sim plotting functions are simply wrappers for ggplot2 code, specific to the output from ss3sim get_results_all() objects. They are designed to quickly explore simulation output, rather than produce publication-level figures. The functions use arguments passed as characters that refer to columns of data. Scalar plots requires a value for x; while, for time-series plots, x = "year" will be necessary.

Note that there are some subtle differences between the functions. Boxplots cannot have a color mapped to them like points or lines, and thus, color is not a valid argument. The time-series point and line plots are grouped internally by 'ID', which is a combination of scenario and iteration and will be automatically added to the data set if not already present.

Output

These functions print the ggplot object, but also return it invisibly for saving or printing again later. For example, you could save the ggplot object and add a custom theme or change an axis label before printing it.

Author(s)

Cole Monnahan

Examples

# Plot scalar values
data("scalar_dat", package = "ss3sim")
re <- calculate_re(scalar_dat)
## Not run: 
plot_boxplot(re,
  x = "E", y = "depletion_re", horiz = "D",
  relative.error = TRUE
)

## End(Not run)
rm(re)
# Merge scalar and time-series values to plot time series with color
# Be patient, the time-series boxplots take some time.
data("ts_dat", package = "ss3sim")
ts_dat[, "model_run"] <- factor(ts_dat[, "model_run"],
  levels = c("om", "em")
)
## Not run: 
plot_boxplot(ts_dat,
  x = "year", y = "SpawnBio",
  horiz = "scenario", vert = "model_run"
)

## End(Not run)

Description

Plot the cumulative mean for a parameter

Usage

plot_cummean(
  data,
  var,
  order_var = "iteration",
  group = NULL,
  use_facet = FALSE
)

Arguments

Value

A list containing the ggplot object and the data used to make it

Examples

data("scalar_dat", package = "ss3sim")
obj <- plot_cummean(scalar_dat[scalar_dat$model_run == "em", ],
  "VonBert_K_Fem_GP_1",
  group = "scenario",
  use_facet = TRUE
)
# obj$plot
# obj$data
rm(obj)
# group can also be left NULL if only plotting a single scenario.
# it is recommended to set use_facet FALSE in this case.
obj2 <- plot_cummean(
  scalar_dat[
    scalar_dat$scenario == unique(scalar_dat$scenario)[1] &
      scalar_dat$model_run == "em",
  ],
  var = "VonBert_K_Fem_GP_1",
  group = NULL,
  use_facet = FALSE
)
# obj2$plot
# obj2$data
rm(obj2)

Description

Plot time-series values as lines

Usage

plot_lines(
  data,
  x = "year",
  y,
  horiz = NULL,
  horiz2 = NULL,
  vert = NULL,
  vert2 = NULL,
  relative.error = FALSE,
  color = NULL,
  axes.free = TRUE,
  print = TRUE
)

Arguments

Details

The ss3sim plotting functions are simply wrappers for ggplot2 code, specific to the output from ss3sim get_results_all() objects. They are designed to quickly explore simulation output, rather than produce publication-level figures. The functions use arguments passed as characters that refer to columns of data. Scalar plots requires a value for x; while, for time-series plots, x = "year" will be necessary.

Note that there are some subtle differences between the functions. Boxplots cannot have a color mapped to them like points or lines, and thus, color is not a valid argument. The time-series point and line plots are grouped internally by 'ID', which is a combination of scenario and iteration and will be automatically added to the data set if not already present.

Output

These functions print the ggplot object, but also return it invisibly for saving or printing again later. For example, you could save the ggplot object and add a custom theme or change an axis label before printing it.

Author(s)

Cole Monnahan

Examples

data("scalar_dat", "ts_dat", package = "ss3sim")
# Merge in max_grad, a performance metric, to use for color
re <- merge(
  by = "ID",
  calculate_re(ts_dat, add = FALSE),
  calculate_re(scalar_dat, add = FALSE)[, c("ID", "max_grad")]
)
## Not run: 
plot_lines(re,
  y = "SpawnBio_re", horiz = "D", vert = "E",
  relative.error = TRUE, color = "max_grad"
)

## End(Not run)

Description

Generate a scatterplot using ggplot2::ggplot to visualize the relationship between two continuous variables.

Usage

plot_points(
  data,
  x,
  y,
  horiz = NULL,
  horiz2 = NULL,
  vert = NULL,
  vert2 = NULL,
  jitter.height = 0,
  jitter.width = 0,
  color = NULL,
  relative.error = FALSE,
  axes.free = TRUE,
  print = TRUE
)

Arguments

Details

Points are placed on the figure using the width setting in ggplot2::position_jitter() that defaults to 40% resolution of the data, meaning that the jitter values will occupy 80% of the implied bins. The previous information was found in the documentation for ggplot2::position_jitter().

Values of NA are removed prior to plotting such that the typical error message from ggplot2 is not printed to the screen.

The ss3sim plotting functions are simply wrappers for ggplot2 code, specific to the output from ss3sim get_results_all() objects. They are designed to quickly explore simulation output, rather than produce publication-level figures. The functions use arguments passed as characters that refer to columns of data. Scalar plots requires a value for x; while, for time-series plots, x = "year" will be necessary.

Note that there are some subtle differences between the functions. Boxplots cannot have a color mapped to them like points or lines, and thus, color is not a valid argument. The time-series point and line plots are grouped internally by 'ID', which is a combination of scenario and iteration and will be automatically added to the data set if not already present.

Output

These functions print the ggplot object, but also return it invisibly for saving or printing again later. For example, you could save the ggplot object and add a custom theme or change an axis label before printing it.

Author(s)

Cole Monnahan

Examples

# Plot scalar values
data("scalar_dat", package = "ss3sim")
re <- calculate_re(scalar_dat)
## Not run: 
plot_points(re,
  x = "E", y = "depletion_re", horiz = "D",
  color = "max_grad", relative.error = TRUE
)

## End(Not run)
rm(re)
# Merge scalar and time-series values to plot time series with color
data("ts_dat", package = "ss3sim")
re <- merge(
  by = "ID",
  calculate_re(ts_dat, add = FALSE),
  calculate_re(scalar_dat, add = FALSE)[, c("ID", "max_grad")]
)
## Not run: 
plot_points(re,
  x = "year", y = "SpawnBio_re",
  horiz = "scenario", color = "max_grad", relative.error = TRUE
)

## End(Not run)
rm(re)

Description

Use ggplot2::ggplot to plot data from ss3sim simulation.

Usage

plot_ss3sim(
  data,
  x,
  y,
  color = NULL,
  relative.error = FALSE,
  axes.free = TRUE,
  print = TRUE,
  horiz = NULL,
  horiz2 = NULL,
  vert = NULL,
  vert2 = NULL
)

Arguments

Details

The ss3sim plotting functions are simply wrappers for ggplot2 code, specific to the output from ss3sim get_results_all() objects. They are designed to quickly explore simulation output, rather than produce publication-level figures. The functions use arguments passed as characters that refer to columns of data. Scalar plots requires a value for x; while, for time-series plots, x = "year" will be necessary.

Note that there are some subtle differences between the functions. Boxplots cannot have a color mapped to them like points or lines, and thus, color is not a valid argument. The time-series point and line plots are grouped internally by 'ID', which is a combination of scenario and iteration and will be automatically added to the data set if not already present.

Output

These functions print the ggplot object, but also return it invisibly for saving or printing again later. For example, you could save the ggplot object and add a custom theme or change an axis label before printing it.

Author(s)

Cole Monnahan

Description

Runs an operating model over a range of fishing mortality, $F$, levels to determine the F at maximum sustainable yield, $F_{MSY}$.

Usage

profile_fmsy(
  om_in,
  results_out,
  start = 0,
  end = 1.5,
  by_val = 0.01,
  verbose = FALSE
)

Arguments

Details

profile_fmsy() runs the operating model with a constant level of fishing for each year and extracts the expected catch in the terminal year. It is assumed that the model time series is long enough for the population to come to equilibrium, and thus, the catch in the terminal year is equivalent to equilibrium catch.

If the function is run with verbose = TRUE, which is not the default, the coefficient of variations of the catches in the terminal years of the model will be printed to the screen. Here, terminal is defined as half as many years as there are ages in the population dynamics of your model. Thus, if the population plus group starts at age twenty, the standard deviation of the last ten years of catch divided by the mean catch over that same time will be printed to the screen for each model that is ran. For the default cod model provided within the package, the CV is less than 1e-04 for all explored levels of fishing mortality.

Ensure that the argument om_in leads to an operating model that is configured for use within ss3sim. For example, the $F$ type must allow for an input vector of $Fs$ rather than catches, along with other specifications.

Value

A data frame of catch by fishing mortality is returned invisibly and saved to the disk along with a figure, Fmsy.pdf.

Examples

## Not run: 
d <- system.file("extdata", "models", "cod-om", package = "ss3sim")
fmsy.val <- profile_fmsy(
  om_in = d, results_out = "fmsy",
  start = 0.1, end = 0.2, by_val = 0.05
)
# cleanup
unlink("fmsy", recursive = TRUE)

## End(Not run)

Description

Replace items with a zero length in a list with the value supplied in the argument replacement. Useful for scenarios where NULL is sometimes a valid input from legacy code. Nested lists can behave badly when they have a NULL entry, for example when converting to a tibble::tibble(), it will be unnamed.

Usage

replace_x(x, replacement = NA_integer_)

Arguments

Value

The object x is returned with some items replaced. If the input object was of zero length, then the replacement parameter will be returned instead.

Author(s)

Amanda from stack overflow

Examples

employees <- list(
  list(
    id = 1,
    dept = "IT",
    age = 29,
    sportsteam = "softball"
  ),
  list(
    id = 2,
    dept = "IT",
    age = 30,
    sportsteam = NULL
  ),
  list(
    id = 3,
    dept = "IT",
    age = 29,
    sportsteam = "hockey"
  ),
  list(
    id = 4,
    dept = NULL,
    age = 29,
    sportsteam = "softball"
  )
)
# Meat of the example here!
## Not run: 
do.call(rbind, lapply(employees, rbind)) %>%
  data.frame() %>%
  purrr::modify_depth(2, replace_x)

## End(Not run)

Description

This is the main high-level wrapper function used to run a set of ss3sim simulations. The data frame passed to simdf is parsed into a list and used to control ss3sim_base(). Alternatively, you can call ss3sim_base() directly with your own lists.

Usage

run_ss3sim(
  iterations,
  simdf = NULL,
  parallel = FALSE,
  parallel_iterations = FALSE,
  ...
)

Arguments

library(doParallel)
cl <- makeCluster(2)
registerDoParallel(cl)

Details

The operating model folder, which is passed as a file path using simdf[["om_dir"]], should contain the following files:

• forecast.ss,

• yourmodel.ctl,

• yourmodel.dat,

• ss.par, and

• starter.ss. The files should be the formatted versions that are returned from Stock Synthesis after the model is optimized, i.e., .ss_new files. It is important to use these formatted files because many functions in ss3sim and r4ss depend on the location of keywords present in the comments and other standardized formatting. Once you have these files from a successfully optimized model, rename the .ss_new files to match the names listed above, though you can change yourmodel to whatever name is listed for the control and data files in starter.ss. The estimation model folder should also contain these files, except ss.par and yourmodel.dat files, which are unnecessary. See the vignette titled modifying-models for details on modifying an existing Stock Synthesis model to run with ss3sim. Alternatively, consider modifying the built-in model configuration based on north sea cod.

Note that due to the way that Stock Synthesis is being used as an OM, you may see the following error from ADMB in the console: Error – base = 0 in function prevariable& pow(const prevariable& v1, CGNU_DOUBLE u) However, this is not a problem because ADMB is not used to optimize the OM, and thus, the error can safely be ignored.

Value

The output will appear in your current R working directory. Folders will be named based on the "scenario" column of simdf or based on the date-time stamp (i.e., mmddhhmmss) generated automatically at the start of the simulation. The resulting folders will look like the following if you run your simulation at noon on January 01:

• 0101120000/1/om

• 0101120000/1/em

• 0101120000/2/om

• ...

Author(s)

Sean C. Anderson

See Also

ss3sim_base() can be called directly by passing lists to each individual argument rather than using the data-frame approach of run_ss3sim(simdf = ). The lists correspond to each function called by ss3sim_base(). Each element is itself a list of arguments for the given function. Either way allows users to pass arguments to each of the ⁠change_*()⁠ or ⁠sample_*()⁠ functions. Note that if you do not include an argument, then ss3sim_base() will assume the value of that argument is NULL.

Examples

## Not run: 
# A run with deterministic process error for model checking
# by passing user_recdevs to ss3sim_base through run_ss3sim:
recdevs_det <- matrix(0, nrow = 101, ncol = 2)
df <- data.frame(setup_scenarios_defaults(),
  "scenarios" = "determinate"
)
run_ss3sim(
  iterations = 1:2, simdf = df,
  bias_adjust = FALSE, user_recdevs = recdevs_det
)
get_results_all(user_scenarios = "determinate", overwrite = TRUE)
ts <- utils::read.csv("ss3sim_ts.csv")
expect_equivalent(
  unlist(ts$rec_dev[ts$year %in% 1:10 & ts$iteration == 2]),
  recdevs_det[1:10, 2]
)

## End(Not run)

Description

Extract age-composition data from a .ss_new data file and sample the data. It is assumed that the composition data will be expected values as written by Stock Synthesis in the second section of the data file, but one can also sample input data. The resulting age-composition data are assumed to represent observed age composition and will overwrite the age data in dat_list, which is returned invisibly. The data file can also be written to the disk, if a file path is provided to outfile, and used as simulated data by an estimation model.

Usage

sample_agecomp(
  dat_list,
  outfile = NULL,
  fleets,
  Nsamp,
  years,
  cpar = 1,
  ESS = NULL,
  keep_conditional = TRUE,
  ...
)

Arguments

Value

A modified .dat file if !is.null(outfile). A list object containing the modified .dat file is returned invisibly.

Author(s)

Cole Monnahan and Kotaro Ono

See Also

Other sampling functions: clean_data(), sample_calcomp(), sample_catch(), sample_discard(), sample_index(), sample_lcomp(), sample_mlacomp(), sample_wtatage()

Examples

d <- system.file("extdata", package = "ss3sim")
f_in <- file.path(d, "models", "cod-om", "codOM.dat")
dat_list <- r4ss::SS_readdat(f_in, verbose = FALSE)

## Turn off age comps by specifying fleets=NULL
test <- sample_agecomp(dat_list = dat_list, fleets = NULL)

## Generate with a smaller number of fleet taking samples
ex1 <- sample_agecomp(
  dat_list = dat_list, outfile = NULL,
  fleets = 2, Nsamp = list(c(10, 50)), years = list(c(26, 27))
)
NROW(ex1$agecomp) == 2

## Generate with varying Nsamp by year for first fleet
ex2 <- sample_agecomp(
  dat_list = dat_list, outfile = NULL,
  fleets = c(1, 2),
  Nsamp = list(c(rep(50, 5), rep(100, 5)), 50),
  years = list(seq(26, 44, 2), c(26:100))
)

## Run three  cases showing Multinomial, Dirichlet(1), and over-dispersed
## Dirichlet for different levels of sample sizes
op <- graphics::par(mfrow = c(1, 3))
set.seed(1)
true <- prop.table(dat_list$agecomp[
  dat_list$agecomp$FltSvy == 1 & dat_list$agecomp$Yr == 50, -(1:9)
])
cpars <- c(NA, 1, 4)
for (samplesize in c(30, 100, 1000)) {
  if (samplesize > 30) graphics::par(mar = c(5.1, 1, 4.1, 2.1))
  graphics::plot(dat_list$agebin_vector, true,
    type = "b", ylim = c(0, 1),
    col = 4, lwd = 2, xlab = "Age",
    ylab = ifelse(samplesize == 30, "Proportion", ""),
    main = paste("Sample size =", samplesize)
  )
  if (samplesize == 30) {
    graphics::legend("topright",
      lty = 1, col = 1:4, bty = "n",
      legend = c("Multinomial", "Dirichlet(1)", "Dirichlet(4)", "Truth")
    )
  }
  for (i in seq_along(cpars)) {
    ex <- sample_agecomp(
      dat_list = dat_list, outfile = NULL, fleets = 1,
      Nsamp = list(samplesize), years = list(50), cpar = cpars[i]
    )$agecomp
    lines(dat_list$agebin_vector, prop.table(ex[1, -(1:9)]),
      col = i, type = "b"
    )
  }
}
graphics::par(op)

Description

Sample conditional age-at-length (CAAL) data from expected values of length proportions and expected values of age proportions (conditional on length) from the operating model (OM) and writes the samples to file for use by the estimation model (EM).

Usage

sample_calcomp(
  dat_list,
  exp_vals_list,
  outfile = NULL,
  fleets,
  years,
  Nsamp_lengths,
  Nsamp_ages,
  method = "simple_random",
  ESS_lengths = NULL,
  ESS_ages = NULL,
  lcomps_sampled = FALSE,
  ...
)

Arguments

Details

There are many steps needed to sample CAAL data because ages are not independent from lengths. The data is located in the .dat file alongside age compositions. CAAL have the added complexity of one line per length bin. Thus, each row represents the observed age distribution for a length bin conditioned on the fish lengths that were observed in the length compositions. The age distribution will be truncated for older or younger fish. Often, many rows will be empty because no fish of that length bin were observed. These empty rows are not needed in the .dat file.

The sampling process includes the following steps:

1. Lengths are sampled based on the desired number of lengths, $N$. $N$ is the maximum amount that could be aged.

2. Those lengths are binned to create a length distribution, i.e., numbers of fish in each length bin.

3. Ages are sampled from fish that contributed to the length distribution. Several strategies are possible for sampling ages from those fish

4. age all fish,

5. take random subset of fish independent of length bin, or

6. take a fixed number of fish from each length bin.

ss3sim can currently only handle randomly sampling ages from lengthed fish. Future versions could include the last option; please contact the developers if you are interested in helping facilitate this.

Note that the overall total sample size for all CAAL bins is specified by the user for the given fleet and year in Nsamp_ages. These sample sizes and the expected values of age proportions (conditional on length) are used to sample for realistic age proportions. If all fish are aged, then no resampling is performed. If no fish are aged for a row of age proportions in conditional age at length data, then that row is discarded. If all fish are not aged, then a new sample size must be drawn. This new sample size must be less than or equal to the number of fish that were sampled for their length. This new sample size is used to draw ages randomly from the expected values. If we consider all rows for a fleet and year (one for each length bin), then the sum of those will be the sample size for the CAAL data. However, if the CAAL sample size is less than the length sample size, We accomplish this in the code by doing sampling without replacement for vectors of length bins equal to the number of fish in them. This ensures realistic sampling. If the option (3) above were implemented, a different strategy would need to be implemented. For instance, if the user wants 10 fish from each length bin but only 5 fish were observed, what to do? A value of NULL for fleets indicates to delete the CAAL data but not the marginal age data.

When Dirichlet sampling is used for length compositions, the number of fish observed will be real-valued and not whole fish. One cannot simply multiply by the length composition sample size to get whole numbers because they are real and rounding or truncating would be unsatisfactory. Currently, the function simply draws a multinomial sample from the length compositions of specified size (Nsamp). However, this does not guarantee that fewer fish are aged than lengthed. If you are specifying a small number of fish to age relative to length, then this might be alright. However, we discourage the use of Dirichlet length samples when using CAAL data as currently implemented.

Note that this function cannot handle all types of CAAL sampling. This function requires that there be a row of CAAL data for each length data bin (for each year and fleet that sampling is specified to be performed), where Lbin_lo and Lbin_hi are the same value. Note also that this sampling procedure represents simple random sampling for CAAL, where (1) lengths are sampled randomly, (2) fish are lengthed and placed into bins, and (3) a subset of lengthed fish are aged, where a constant proportion from each length bin are selected for aging. This does not represent length stratified sampling where a subset of lengthed fish are aged, and a constant number from each length bin is selected for aging, although these data could also be put into a Stock Synthesis model as CAAL.

Value

A modified .dat file if !is.null(outfile). A list object containing the modified .dat file is returned invisibly.

Author(s)

Cole Monnahan, Kotaro Ono

See Also

Other sampling functions: clean_data(), sample_agecomp(), sample_catch(), sample_discard(), sample_index(), sample_lcomp(), sample_mlacomp(), sample_wtatage()

Description

This function creates a matrix of sampled catches from the expected available catches for all fleets with catches. The input value used for catch_se will be used to resample the catches. There is a bit of a disconnect here because catches are defined by input F values not absolute catches. Let $D_y$ be the discard from the operating model for year y. Then the sampled value is calculated as: $D_y*exp(stats::rnorm(1, 0, sds_obs)-sds_obs^2/2)$. The second term adjusts the random samples so that their expected value is $D_y$, i.e., the log-normal bias correction.

Usage

sample_catch(dat_list, outfile = NULL)

Arguments

Value

A modified .dat file if !is.null(outfile). A list object containing the modified .dat file is returned invisibly.

Author(s)

Kelli F. Johnson

See Also

Other sampling functions: clean_data(), sample_agecomp(), sample_calcomp(), sample_discard(), sample_index(), sample_lcomp(), sample_mlacomp(), sample_wtatage()

Description

Apply the multinomial or Dirichlet distribution to sample composition data, creating a data frame that mimics observed composition data.

Usage

sample_comp(data, Nsamp, fleets, years, ESS = NULL, cpar = 1, ...)

Arguments

Details

Sample size, i.e., 'Nsamp', is used as a measure of precision, where higher sample sizes lead to simulated samples that more accurately represent the truth provided in data.

Value

A data frame of observed composition data.

Author(s)

Kelli F. Johnson

Description

This function creates an index of discards sampled from the expected available discards for specified fleets in specified years. Let $D_y$ be the discard from the operating model for year y. Then the sampled value is calculated as: $D_y*exp(stats::rnorm(1, 0, sds_obs)-sds_obs^2/2)$. The second term adjusts the random samples so that their expected value is $D_y$, i.e., the log-normal bias correction.

Usage

sample_discard(
  dat_list,
  outfile = NULL,
  fleets,
  years,
  sds_obs,
  seas = list(1)
)

Arguments

Value

A modified .dat file if !is.null(outfile). A list object containing the modified .dat file is returned invisibly.

Author(s)

Kelli F. Johnson

See Also

Other sampling functions: clean_data(), sample_agecomp(), sample_calcomp(), sample_catch(), sample_index(), sample_lcomp(), sample_mlacomp(), sample_wtatage()

Description

Sample with a Dirichlet-Multinomial distribution

Usage

sample_dm(data, n, par)

Arguments

Value

A data frame with one row because right now the input data should only be a single row of data.

Description

Create new catch-per-unit-effort (CPUE)/indices of abundance that are based on the numbers in a data file. Typically the data file will be filled with expected values rather than observed data but it does not have to be. Sampling can only occur on fleets, years, and seasons that have current observations. If rows of information are not sampled from, then they are removed. So, you can take away rows of data but you cannot add them with this function.

Usage

sample_index(
  dat_list,
  outfile = lifecycle::deprecated(),
  fleets,
  years,
  sds_obs = list(0.01),
  sds_out,
  seas = list(1)
)

Arguments

Details

Limitations to the functionality of this function are as follows:

• you can only generate observations from rows of data that are present, e.g., you cannot make a new observation for a year that is not present in the passed data file;

• no warning will be given if some of the desired year, seas, fleet combinations are available but not all, instead just the combinations that are available will be returned in the data list object; and

• sampling uses a log-normal distribution when the log-normal distribution is specified in CPUEinfo[["Errtype"]] and a normal distribution for all other error types, see below for details on the log-normal sampling.

Samples are generated using the following equation when the log-normal distribution is specified:

$B_y*exp(stats::rnorm(1, 0, sds_obs)-sds_obs^2/2)$

, where $B_y$ is the expected biomass in year y and $sds_obs$ is the standard deviation of the normally distributed biomass or the standard error of the $log_e(B_y)$. For the error term, this is the same parameterization that is used in Stock Synthesis. More details can be found in the section on indices in the Stock Synthesis manual The second term in the equation adjusts the random samples so their expected value is $B_y$, i.e., the log-normal bias correction.

If you only know the coefficient of variation ($CV$), then the input error can be approximated using $\sqrt{log_e(1+CV^{2})}$. Where, $CV$ is assumed to be constant with mean changes in biomass. The log-normal distribution can be approximated by a proportional distribution or normal distribution only when the variance is low, i.e., $CV < 0.50$ or log standard deviation of 0.22.

Value

A Stock Synthesis data file list object is returned. The object will be a modified version of dat_list.

Author(s)

Cole Monnahan, Kotaro Ono

See Also

Other sampling functions: clean_data(), sample_agecomp(), sample_calcomp(), sample_catch(), sample_discard(), sample_lcomp(), sample_mlacomp(), sample_wtatage()

Examples

# Add a list from [r4ss::SS_readdat()] to your workspace, this is example
# data that is saved in the ss3sim package.
# Index data are saved in `dat_list[["CPUE"]]`
dat_list <- r4ss::SS_readdat(
  file = file.path(
    system.file("extdata", "example-om", package = "ss3sim"),
    "ss3_expected_values.dat"
  ),
  verbose = FALSE
)


# Sample from each available year from fleet 2 with an increasing trend in
# the observation error, i.e., the most recent year has the highest
# likelihood to be the furthest from the input data
ex1 <- sample_index(
  dat_list,
  outfile = NULL,
  fleets = 2,
  seas = list(
    dat_list[["CPUE"]][dat_list[["CPUE"]][, "index"] == 2, "seas"]
  ),
  years = list(dat_list[["CPUE"]][["year"]]),
  sds_obs = list(
    seq(0.001, 0.1, length.out = length(dat_list[["CPUE"]][["year"]]))
  )
)


## Not run: 
# Sample from less years, note that sampling from more years than what is
# present in the data will not work
ex2 <- sample_index(dat_list,
  outfile = NULL,
  fleets = 2,
  seas = list(unique(
    dat_list[["CPUE"]][dat_list[["CPUE"]][, "index"] == 2, "seas"]
  )),
  years = list(dat_list[["CPUE"]][["year"]][-c(1:2)]),
  sds_obs = list(0.001)
)

# sd in the returned file can be different than what is used to sample, this
# is helpful when you want to test what would happen if the estimation method
# was improperly specified
ex3 <- sample_index(
  dat_list = dat_list,
  fleets = 2,
  seas = list(unique(
    dat_list[["CPUE"]][dat_list[["CPUE"]][, "index"] == 2, "seas"]
  )),
  years = list(dat_list[["CPUE"]][["year"]]),
  sds_obs = list(0.01),
  sds_out = list(0.20)
)
ex3[["CPUE"]][["se_log"]]

## End(Not run)
# Sample from two fleets after adding fake CPUE data for fleet 1
dat_list2 <- dat_list
dat_list2[["CPUE"]] <- rbind(
  dat_list[["CPUE"]],
  dat_list[["CPUE"]] |>
    dplyr::mutate(index = 1, seas = 1)
)
dat_list2[["N_cpue"]] <- NROW(dat_list2[["CPUE"]])
ex4 <- sample_index(
  dat_list = dat_list2,
  fleets = 1:2,
  seas = list(1, 7),
  # Subset two years from each fleet
  years = list(c(76, 78), c(80, 82)),
  # Use the same sd values for both fleets
  sds_obs = list(0.01),
  sds_out = list(0.20)
)

Description

Extract length-composition data from a .ss_new data file and sample the data. It is assumed that the composition data will be expected values as written by Stock Synthesis in the second section of the data file, but one can also sample input data. The resulting length-composition data are assumed to represent observed length composition and will overwrite the length data in dat_list, which is returned invisibly. The data file can also be written to the disk, if a file path is provided to outfile, and used as simulated data by an estimation model.

Usage

sample_lcomp(
  dat_list,
  outfile = NULL,
  fleets,
  Nsamp,
  years,
  cpar = 1,
  ESS = NULL,
  ...
)

Arguments

Value

A modified .dat file if !is.null(outfile). A list object containing the modified .dat file is returned invisibly.

Author(s)

Cole Monnahan and Kotaro Ono

See Also

sample_agecomp() for more examples.

Other sampling functions: clean_data(), sample_agecomp(), sample_calcomp(), sample_catch(), sample_discard(), sample_index(), sample_mlacomp(), sample_wtatage()

Examples

dat_list <- r4ss::SS_readdat(
  verbose = FALSE,
  file = system.file(file.path("extdata", "models", "cod-om", "codOM.dat"),
    package = "ss3sim"
  )
)
## Generate with constant sample size across years
ex1 <- sample_lcomp(
  dat_list = dat_list, outfile = NULL,
  fleets = 1:2, Nsamp = list(100, 50),
  years = list(seq(26, 100, by = 2), 80:100)
)

Description

Sample a standard normal in log-space and apply the error to observations.

Usage

sample_lognormal(obs, sd)

Arguments

Details

Newly sampled values are calculated $obs*exp(stats::rnorm(1, 0, sd)-sd^2/2)$. The second term adjusts the random samples so that their expected value is obs, i.e., the log-normal bias correction.