Skip to content

Commit

Permalink
WIP: Base resampling model
Browse files Browse the repository at this point in the history 
initial attempt to implement a base model that is used to sample historic expected values for OM and bootstrap values for EM
  • Loading branch information
@nathanvaughan-NOAA
nathanvaughan-NOAA committed Jan 5, 2025
1 parent 3d9eefa commit 1c50dd0
Show file tree
Hide file tree
Showing 3 changed files with 215 additions and 9 deletions.
110 changes: 110 additions & 0 deletions R/initOM.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -505,6 +505,116 @@ run_OM <- function(OM_dir,
return(dat)
}


#' Initial run of the Base model to get expected and bootstrap values
#'
#' This function is used to initialize the Base model and get expected values
#' for the OM and bootstraps for the EM.
#' @template OM_dir
#' @template verbose
#' @param debug_par_run If set to TRUE, and the run fails, a new folder called
#' error_check will be created, and the model will be run from control start
#' values instead of ss.par. The 2 par files are then compared to help debug
#' the issue with the model run. Defaults to TRUE.
#' @template seed
#' @author Kathryn Doering
#' @importFrom r4ss SS_readdat SS_readstarter SS_writestarter
run_Base <- function(OM_dir,
verbose = FALSE,
debug_par_run = TRUE,
seed = NULL) {


if (is.null(seed)) {
seed <- stats::runif(1, 1, 9999999)
}

start <- r4ss::SS_readstarter(file.path(OM_dir, "starter.ss"),
verbose = FALSE
)

start[["init_values_src"]] <- 1
start[["N_bootstraps"]] <- 2
start[["seed"]] <- seed

r4ss::SS_writestarter(start,
dir = OM_dir, verbose = FALSE, overwrite = TRUE
)

# run SS and get the data set
run_ss_model(OM_dir, "-maxfn 0 -phase 50 -nohess",
verbose = verbose,
debug_par_run = debug_par_run
)

dat <- r4ss::SS_readdat(file.path(OM_dir, "data.ss_new"),
section = 2,
verbose = FALSE
)

return(dat)
}

#' Refit OM model using expected values from original OM
#'
#' This function is used to initialize the Base model and get expected values
#' for the OM and bootstraps for the EM.
#' @template OM_dir
#' @param base_dat Expected Value data list from base model run.
#' @template verbose
#' @param debug_par_run If set to TRUE, and the run fails, a new folder called
#' error_check will be created, and the model will be run from control start
#' values instead of ss.par. The 2 par files are then compared to help debug
#' the issue with the model run. Defaults to TRUE.
#' @template seed
#' @author Kathryn Doering
#' @importFrom r4ss SS_readdat SS_readstarter SS_writestarter
refit_OM <- function(OM_dir,
base_dat,
verbose = FALSE,
debug_par_run = TRUE,
seed = NULL) {

if (is.null(seed)) {
seed <- stats::runif(1, 1, 9999999)
}

start <- r4ss::SS_readstarter(file.path(OM_dir, "starter.ss"),
verbose = FALSE
)
start[["init_values_src"]] <- 1
start[["detailed_age_structure"]] <- 1
start[["last_estimation_phase"]] <- 10
start[["depl_basis"]] <- 1
start[["depl_denom_frac"]] <- 1
start[["SPR_basis"]] <- 4
start[["F_report_units"]] <- 1
start[["F_report_basis"]] <- 0
start[["F_age_range"]] <- NULL
start[["ALK_tolerance"]] <- 0
start[["minyr_sdreport"]] <- -1
start[["maxyr_sdreport"]] <- -2
start[["N_bootstraps"]] <- 3
start[["seed"]] <- seed

r4ss::SS_writestarter(start,
dir = OM_dir, verbose = FALSE, overwrite = TRUE
)

r4ss::SS_writedat(base_dat, file.path(OM_dir, start[["datfile"]]),
overwrite = TRUE,
verbose = FALSE
)

# run SS and get the data set
run_ss_model(OM_dir, "-nohess",
verbose = verbose,
debug_par_run = debug_par_run
)

return(NULL)
}

#' Get the sampling scheme in a data file.
#'
#' Determine what the default sampling scheme is for a given data file.
Expand Down
43 changes: 41 additions & 2 deletions R/runSSMSE.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -556,23 +556,41 @@ run_SSMSE_iter <- function(out_dir = NULL,
OM_in_dir <- out_loc[["OM_in_dir"]]
EM_in_dir <- out_loc[["EM_in_dir"]]
EM_out_dir <- out_loc[["EM_out_dir"]]
Base_out_dir <- out_loc[["Base_out_dir"]]
Base_in_dir <- out_loc[["Base_in_dir"]]
if (!is.null(EM_out_dir)) {
EM_out_dir_basename <- strsplit(EM_out_dir, "_init$")[[1]][1]
} else {
EM_out_dir_basename <- NULL
}
copy_model_files(
OM_in_dir = OM_in_dir, OM_out_dir = OM_out_dir,
EM_in_dir = EM_in_dir, EM_out_dir = EM_out_dir
EM_in_dir = EM_in_dir, EM_out_dir = EM_out_dir,
Base_in_dir = Base_in_dir, Base_out_dir= Base_out_dir
)

# clean model files ----
# want to do this as soon as possible to "fail fast"
# for now, this just gets rid of -year value observations, but could do
# other things.
clean_init_mod_files(
OM_out_dir = OM_out_dir, EM_out_dir = EM_out_dir,
Base_out_dir= Base_out_dir,
overwrite = TRUE
)
# Add code to create a base model that will be used to sample expected values
# to inform an OM and bootstrapped values to inform the EM's
Base_dat <- run_Base(
OM_dir = Base_out_dir, verbose = verbose,
seed = (iter_seed[["iter"]][1] + 12345)
)

browser()

refit_OM(
OM_dir = OM_out_dir, base_dat = Base_dat,
verbose = verbose, seed = (iter_seed[["iter"]][1] + 12345)
)

# convert sample_struct names ----
# get the full sampling structure for components that the user didnt specify.
Expand All @@ -588,6 +606,22 @@ run_SSMSE_iter <- function(out_dir = NULL,
if(!is.null(sample_struct_hist)){
sample_struct_hist <- convert_to_r4ss_names(sample_struct_hist)
}
# Convert the user input parameter modifications into vectors of annual additive deviations
future_base_dat <- convert_future_om_list_to_devs_df(future_om_list = future_om_list, scen_name = scen_name, niter = niter, om_mod_path = Base_out_dir, nyrs = nyrs, global_seed = (iter_seed[["iter"]][1] + 1234))

# MSE first iteration ----
# turn the stock assessment model into an OM
init_base_mod <- create_OM(
OM_out_dir = Base_out_dir, overwrite = TRUE,
sample_struct_hist = sample_struct_hist,
sample_struct = sample_struct,
verbose = verbose, writedat = TRUE, nyrs = nyrs,
nyrs_assess = nyrs_assess, nscen = nscen,
scen_name = scen_name, niter = niter,
future_om_dat = future_base_dat,
seed = (iter_seed[["iter"]][1] + 1234)
)

# Convert the user input parameter modifications into vectors of annual additive deviations
future_om_dat <- convert_future_om_list_to_devs_df(future_om_list = future_om_list, scen_name = scen_name, niter = niter, om_mod_path = OM_out_dir, nyrs = nyrs, global_seed = (iter_seed[["iter"]][1] + 1234))

Expand All @@ -606,8 +640,13 @@ run_SSMSE_iter <- function(out_dir = NULL,
impl_error <- init_mod[["impl_error"]]
# Complete the OM run so it can be use for expect values or bootstrap
if (use_SS_boot == TRUE) {
# OM_dat <- run_OM(
# OM_dir = OM_out_dir, boot = use_SS_boot, nboot = 1,
# sample_catch = sample_catch,
# verbose = verbose, init_run = TRUE, seed = (iter_seed[["iter"]][1] + 12345)
# )
OM_dat <- run_OM(
OM_dir = OM_out_dir, boot = use_SS_boot, nboot = 1,
OM_dir = Base_out_dir, boot = use_SS_boot, nboot = 1,
sample_catch = sample_catch,
verbose = verbose, init_run = TRUE, seed = (iter_seed[["iter"]][1] + 12345)
)
Expand Down
71 changes: 64 additions & 7 deletions R/utils.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -140,11 +140,12 @@ create_scen_list <- function(scen_name_vec,
#'
#' @template OM_out_dir
#' @template EM_out_dir
#' @param Base_out_dir The directory where the base model files will be written.
#' @template MS
#' @template overwrite
#' @importFrom r4ss SS_readstarter SS_readdat SS_writedat
clean_init_mod_files <- function(OM_out_dir, EM_out_dir = NULL, MS = "EM",
overwrite = FALSE) {
clean_init_mod_files <- function(OM_out_dir, EM_out_dir = NULL, Base_out_dir = NULL,
MS = "EM", overwrite = FALSE) {
# read in files
OM_start <- SS_readstarter(file.path(OM_out_dir, "starter.ss"),
verbose = FALSE
Expand All @@ -157,6 +158,17 @@ clean_init_mod_files <- function(OM_out_dir, EM_out_dir = NULL, MS = "EM",
datlist = OM_dat,
verbose = FALSE
)
Base_start <- SS_readstarter(file.path(Base_out_dir, "starter.ss"),
verbose = FALSE
)
Base_dat <- SS_readdat(file.path(Base_out_dir, Base_start[["datfile"]]),
verbose = FALSE,
section = 1
)
Base_ctl <- SS_readctl(file.path(Base_out_dir, Base_start[["ctlfile"]]),
datlist = Base_dat,
verbose = FALSE
)
if (!is.null(EM_out_dir)) {
EM_start <- SS_readstarter(file.path(EM_out_dir, "starter.ss"),
verbose = FALSE
Expand All @@ -178,9 +190,6 @@ clean_init_mod_files <- function(OM_out_dir, EM_out_dir = NULL, MS = "EM",
styr <- OM_dat[["styr"]]
endyr <- OM_dat[["endyr"]]




# get years in range function
get_yrs_in_range <- function(list_name, dat, styr, endyr) {
df <- dat[[list_name]]
Expand Down Expand Up @@ -513,11 +522,17 @@ create_out_dirs <- function(out_dir, niter, OM_name, OM_in_dir,
# figure out the OM_name and create the directory within the one just created
if (!is.null(OM_name)) {
OM_folder_name <- paste0(OM_name, "_OM")
Base_folder_name <- paste0(OM_name, "_OM_Base")
} else {
OM_folder_name <- paste0(basename(OM_in_dir), "_OM")
Base_folder_name <- paste0(basename(OM_in_dir), "_OM_Base")
}
OM_out_dir <- file.path(out_dir, OM_folder_name)
dir.create(OM_out_dir, showWarnings = FALSE)

Base_out_dir <- file.path(out_dir, Base_folder_name)
dir.create(Base_out_dir, showWarnings = FALSE)

# Add the EM dir, if necessary
if (is.null(EM_name) & !is.null(EM_in_dir)) EM_name <- basename(EM_in_dir)
if (!is.null(EM_name) & is.null(EM_in_dir)) {
Expand All @@ -542,7 +557,8 @@ create_out_dirs <- function(out_dir, niter, OM_name, OM_in_dir,
}
OM_mod_loc <- list(
OM_in_dir = OM_in_dir, OM_out_dir = OM_out_dir,
EM_in_dir = EM_in_dir, EM_out_dir = EM_out_dir
EM_in_dir = EM_in_dir, EM_out_dir = EM_out_dir,
Base_in_dir = OM_in_dir, Base_out_dir = Base_out_dir
)
}

Expand Down Expand Up @@ -587,13 +603,17 @@ locate_in_dirs <- function(OM_name = NULL, OM_in_dir = NULL) {
#' @template OM_EM_in_dir
#' @template OM_out_dir
#' @template EM_out_dir
#' @param Base_in_dir Relative or absolute path to the input base model
#' @param Base_out_dir Relative or absolute path to the base model outputs
#' @template verbose
#' @return TRUE, if copying is successful
#'
copy_model_files <- function(OM_in_dir = NULL,
OM_out_dir = NULL,
EM_in_dir = NULL,
EM_out_dir = NULL,
Base_in_dir = NULL,
Base_out_dir= NULL,
verbose = FALSE) {
# checks
if (!is.null(OM_in_dir)) {
Expand Down Expand Up @@ -656,7 +676,44 @@ copy_model_files <- function(OM_in_dir = NULL,
} else {
success_EM <- TRUE
}
invisible(c(success_OM = success_OM, success_EM = success_EM))

if (!is.null(Base_in_dir)) {
if (!all(c(
"control.ss_new", "data.ss_new", "starter.ss_new",
"forecast.ss_new", "ss.par"
) %in% list.files(Base_in_dir))) {
stop(
".ss_new files not found in the original Base directory ",
Base_in_dir, ". Please run the model to make the .ss_new files available."
)
}
}
# copy over Base ----
if (!is.null(Base_in_dir) & !is.null(Base_out_dir)) {
if (verbose == TRUE) {
message(
"Copying over .ss_new model files in ", Base_in_dir,
" to ", Base_out_dir, "."
)
}
success_Base <- r4ss::copy_SS_inputs(
dir.old = Base_in_dir,
dir.new = Base_out_dir,
overwrite = FALSE,
use_ss_new = TRUE, # will rename the ss new files, also.
copy_par = TRUE,
verbose = FALSE
)
if (success_Base == FALSE) {
stop(
"Problem copying SS Base .ss_new files from ", Base_in_dir, " to ",
Base_out_dir, "."
)
}
} else {
success_Base <- TRUE
}
invisible(c(success_OM = success_OM, success_EM = success_EM, success_Base = success_Base))
}

#' function that creates a combined column to the list_item of interest
Expand Down

0 comments on commit 1c50dd0

Please sign in to comment.