actually commit select_n_knots

e-perl-NOAA · e-perl-NOAA · commit 6c6236a75882 · 2026-02-05T16:14:03.000-05:00
diff --git a/R/clean_and_resample.R b/R/clean_and_resample.R
@@ -27,7 +27,8 @@
 #' @param dir_out A character string specifying the directory for output files.
 #' @param test Logical. TRUE/FALSE. If TRUE, will only run first two resampling
 #' tests.
-#' @param n_knots Numeric. Default  = 500.
+#' @param n_knots Numeric. This will be calculated for you based on sample size 
+#' if there is no input.
 #' @param model_type String. Default = "wrapper_sdmtmb", but can be any preset
 #' wrapper_*() function or a premade home built function.
 #' @param bio A data frame containing the biological data, if applicable. NULL
diff --git a/R/resample_tests.R b/R/resample_tests.R
@@ -37,8 +37,43 @@
 #' }
 #' @examples
 #' \dontrun{
-#' resample_tests() # TO DO: NEED EXAMPLE OF HOW TO USE
+#'   catch <- surveyresamplr::noaa_nwfsc_catch
+#'   grid_yrs <- sdmTMB::replicate_df(
+#'    dat = surveyresamplr::noaa_nwfsc_catch, time_name = "year",
+#'    time_values = unique(catch$year)
+#'   )
+#'  spp_list <- data.frame(
+#'    srvy = "CA",
+#'    common_name = "arrowtooth flounder",
+#'    filter_lat_gt = 34,
+#'    filter_lat_lt = NA,
+#'    filter_depth = NA,
+#'    model_fn = "total_catch_wt_kg ~ 0 + factor(year) + pass",
+#'    model_family = "delta_gamma",
+#'    model_anisotropy = TRUE,
+#'    model_spatiotemporal = "iid, iid"
+#'  )
+#' 
+#'  spp_dfs <- cleanup_by_species(
+#'    spp_info = spp_list,
+#'    catch,
+#'    seq_from = 0.1,
+#'    seq_to = 1,
+#'    seq_by = 0.1,
+#'    tot_dataframes = 91,
+#'    replicate_num = 10,
+#'  )
+#' 
+#'  resample_tests(
+#'    spp_dfs = spp_dfs,
+#'    spp_info = spp_list,
+#'    grid_yrs = grid_yrs,
+#'    dir_out = dir_out,
+#'    test = TRUE,
+#'    model_type = "wrapper_sdmtmb"
+#'  )
 #' }
+#' 
 resample_tests <- function(spp_dfs, spp_info, grid_yrs, dir_out, test = FALSE,
                            parallel = FALSE, n_knots = 500,
                            model_type = "wrapper_sdmtmb") {
@@ -136,14 +171,9 @@ resample_tests <- function(spp_dfs, spp_info, grid_yrs, dir_out, test = FALSE,
   if (parallel) {
     n_cores <- future::availableCores()
     n_workers <- round(n_cores * 0.75)
-    if (Sys.info()['sysname'] == 'Windows') {
-      future::plan(future::multisession, workers = n_workers)
-      message("...Running in parallel with multisession")
-    } else {
-      # Use multicore for Linux/macOS for better performance
-      future::plan(future::multicore, workers = n_workers)
-      message("...Running in parallel with multicore")
-    }
+
+    future::plan(future::multisession, workers = n_workers)
+
     results <- furrr::future_map(
       seq_along(spp_files),
       run_parallel_models,
diff --git a/tests/testthat/test-clean_and_resample.R b/tests/testthat/test-clean_and_resample.R
@@ -0,0 +1,98 @@
+# Test for clean_and_resample function in surveyresamplr
+
+test_that("clean_and_resample runs without error and processes realistic bio and catch data", {
+  # Minimal spp_info row
+  spp_info <- data.frame(
+    srvy = "CA",
+    common_name = "arrowtooth flounder",
+    filter_lat_gt = 34,
+    filter_lat_lt = NA,
+    filter_depth = NA,
+    model_fn = "total_catch_wt_kg ~ 0 + factor(year) + pass",
+    model_family = "delta_gamma",
+    model_anisotropy = TRUE,
+    model_spatiotemporal = "iid, iid",
+    stringsAsFactors = FALSE
+  )
+  
+  # Usable catch data matching noaa_nwfsc_catch structure
+  catch <- data.frame(
+    trawlid = 1:3,
+    common_name = rep("arrowtooth flounder", 3),
+    longitude_dd = c(-124.8186, -124.8533, -125.3244),
+    latitude_dd = c(46.71917, 47.64056, 48.32917),
+    year = 2001:2003,
+    pass = c(1, 2, 1),
+    area_swept_ha = c(1.2, 1.5, 1.3),
+    total_catch_numbers = c(16, 22, 50),
+    total_catch_wt_kg = c(28, 40, 60),
+    depth_m = c(270, 130, 115),
+    srvy = rep("CA", 3),
+    stringsAsFactors = FALSE
+  )
+  
+  # grid_yrs with required columns
+  grid_yrs <- data.frame(
+    longitude_dd = rep(c(-124.8186, -124.8533, -125.3244), each = 3),
+    latitude_dd  = rep(c(46.71917, 47.64056, 48.32917), times = 3),
+    pass         = rep(c(1, 2, 1), times = 3),
+    depth_m      = rep(c(270, 130, 115), times = 3),
+    area_km2     = runif(9, 2, 3), # random area between 2 and 3 km2
+    stringsAsFactors = FALSE
+  )
+
+  # Usable bio data matching noaa_nwfsc_bio structure
+  bio <- data.frame(
+    trawlid = 1:3,
+    common_name = rep("arrowtooth flounder", 3),
+    longitude_dd = c(-124.9829, -124.8533, -125.0133),
+    latitude_dd = c(47.01417, 47.72861, 47.84111),
+    year = 2001:2003,
+    pass = c(1, 2, 1),
+    sex = c("F", "M", "F"),
+    length_cm = c(56, 35, 31),
+    age = c(11, 3, 4),
+    depth_m = c(170, 100, 120),
+    project = rep("NWFSC", 3),
+    srvy = rep("CA", 3),
+    stringsAsFactors = FALSE
+  )
+
+  # Temporary output directory
+  dir_out <- withr::local_tempdir()
+
+  # Mock dependencies
+  stub(clean_and_resample, "cleanup_by_species", function(...) {
+    list(
+      data.frame(source = "rep1", trawlid = 1:3),
+      data.frame(source = "rep2", trawlid = 4:6)
+    )
+  })
+  stub(clean_and_resample, "resample_tests", function(...) {
+    TRUE
+  })
+
+  # Run function and check that it doesn't error
+  expect_no_error(
+    clean_and_resample(
+      spp_info = spp_info,
+      catch = catch,
+      seq_from = 0.1,
+      seq_to = 1,
+      seq_by = 0.1,
+      tot_dataframes = 3,
+      replicate_num = 2,
+      grid_yrs = grid_yrs,
+      dir_out = dir_out,
+      bio = bio
+    )
+  )
+
+  # Optionally, check that bio.csv was written (if bio is provided)
+  bio_file <- file.path(dir_out, "CA_arrowtooth_flounder", "bio.csv")
+  expect_true(file.exists(bio_file))
+
+  # Check that the written bio.csv contains expected column names
+  out_bio <- read.csv(bio_file)
+  expect_true(all(c("trawlid", "common_name", "latitude_dd", "depth_m", "source") %in% names(out_bio)))
+})
diff --git a/tests/testthat/test-cleanup_by_species.R b/tests/testthat/test-cleanup_by_species.R
@@ -0,0 +1,59 @@
+# Test for cleanup_by_species function in surveyresamplr
+
+test_that("cleanup_by_species returns expected structure and filters correctly", {
+  # Minimal spp_info row (arrowtooth flounder, filter latitude > 34)
+  spp_info <- data.frame(
+    srvy = "CA",
+    common_name = "arrowtooth flounder",
+    filter_lat_gt = 34,
+    filter_lat_lt = NA,
+    filter_depth = NA,
+    model_fn = "total_catch_wt_kg ~ 0 + factor(year) + pass",
+    model_family = "delta_gamma",
+    model_anisotropy = TRUE,
+    model_spatiotemporal = "iid, iid",
+    stringsAsFactors = FALSE
+  )
+  
+  # Minimal catch data matching required structure
+  catch <- data.frame(
+    srvy = rep("CA", 5),
+    trawlid = 1:5,
+    common_name = c("arrowtooth flounder", "arrowtooth flounder", "arrowtooth flounder", "other fish", "arrowtooth flounder"),
+    total_catch_numbers = c(100, 120, 110, 50, 130),
+    total_catch_wt_kg = c(10, 20, 30, 5, 40),
+    latitude_dd = c(35, 36, 33, 37, 38),
+    longitude_dd = c(-124.8, -124.9, -124.7, -124.8, -125.3),
+    year = c(2001, 2002, 2001, 2002, 2001),
+    pass = c(1, 2, 1, 2, 1),
+    depth_m = c(270, 130, 115, 200, 150),
+    area_swept_ha = c(1.2, 1.5, 1.3, 1.4, 1.6),
+    stringsAsFactors = FALSE
+  )
+
+  # Run function with reduced effort (for test speed)
+  result <- cleanup_by_species(
+    catch = catch,
+    spp_info = spp_info,
+    seq_from = 0.1,
+    seq_to = 0.3,
+    seq_by = 0.1,
+    tot_dataframes = 11,
+    replicate_num = 4
+  )
+  
+  # Check output is a named list of data frames
+  expect_type(result, "list")
+  expect_true(all(vapply(result, is.data.frame, logical(1))))
+  # Check output data frames have at least the columns expected from catch
+  expect_true(all(c("trawlid", "common_name", "latitude_dd", "depth_m") %in% names(result[[1]])))
+  
+  # Check filter works (should only include common_name == "arrowtooth flounder" and latitude_dd > 34)
+  all_lat <- unlist(lapply(result, function(df) df$latitude_dd))
+  all_names <- unlist(lapply(result, function(df) df$common_name))
+  expect_true(all(all_lat > 34))
+  expect_true(all(all_names == "arrowtooth flounder"))
+  
+  # Check that list names are present
+  expect_true(!is.null(names(result)))
+})
diff --git a/tests/testthat/test-resample_tests.R b/tests/testthat/test-resample_tests.R
@@ -0,0 +1,105 @@
+# Mock species SDM wrapper function for testing
+wrapper_sdmtmb <- function(x, y, z, dir_spp, spp_info, n_knots) {
+  # Return a list with mock components
+  list(
+    fit = list(), # Empty object, to be handled by sdmTMB::tidy and sanity
+    index = data.frame(year = 1:3, index = c(1,2,3))
+  )
+}
+
+# Mock sdmTMB::tidy and sdmTMB::sanity
+sdmTMB::tidy <- function(fit, conf.int = FALSE, effects = NULL) {
+  # Return a simple tidy data frame
+  data.frame(term = "intercept", estimate = 1.23)
+}
+sdmTMB::sanity <- function(fit) {
+  # Return a simple sanity data frame
+  data.frame(sane = TRUE)
+}
+
+# Setup test data
+result <- tibble::tibble(
+  source = c("0.1_2", "0.1_4", "0.1_4", "0.2_2", "0.2_3", "0.2_4", "0.2_4", "0.3_2", "0.3_3"),
+  trawlid = c(5, 1, 2, 5, 5, 1, 2, 5, 5),
+  RandomAssignment = rep(1, 9),
+  srvy = rep("CA", 9),
+  common_name = rep("arrowtooth flounder", 9),
+  total_catch_numbers = c(130, 100, 120, 130, 130, 100, 120, 130, 130),
+  total_catch_wt_kg = c(40, 10, 20, 40, 40, 10, 20, 40, 40),
+  latitude_dd = c(38, 35, 36, 38, 38, 35, 36, 38, 38),
+  longitude_dd = c(-125.3, -124.8, -124.9, -125.3, -125.3, -124.8, -124.9, -125.3, -125.3),
+  year = c(2001, 2001, 2002, 2001, 2001, 2001, 2002, 2001, 2001),
+  pass = c(1, 1, 2, 1, 1, 1, 2, 1, 1),
+  depth_m = c(150, 270, 130, 150, 150, 270, 130, 150, 150),
+  area_swept_ha = c(1.6, 1.2, 1.5, 1.6, 1.6, 1.2, 1.5, 1.6, 1.6)
+)
+spp_dfs <- split(result, result$source)
+
+
+spp_info <- data.frame(
+    srvy = "CA",
+    common_name = "arrowtooth flounder",
+    filter_lat_gt = 34,
+    filter_lat_lt = NA,
+    filter_depth = NA,
+    model_fn = "total_catch_wt_kg ~ 0 + factor(year) + pass",
+    model_family = "delta_gamma",
+    model_anisotropy = TRUE,
+    model_spatiotemporal = "iid, iid",
+    stringsAsFactors = FALSE
+  )
+
+grid_yrs <- data.frame(
+  longitude_dd = c(-124.81, -124.85, -125.32, -125.02, -124.55, -123.99, -125.76, -125.46, -124.36, -124.08),
+  latitude_dd  = c(46.85, 47.60, 48.25, 47.81, 46.42, 45.59, 48.08, 47.94, 43.75, 44.67),
+  pass         = rep(0, 10),
+  depth_m      = c(-159, -112, -97, -30, -191, -35, -655, -335, -110, -11),
+  area_km2     = rep(0, 10),
+  srvy         = rep("CA", 10)
+)
+
+dir_out <- tempdir()
+
+test_that("resample_tests runs and creates output files", {
+  # Remove any prior test directory
+  test_dir <- paste0(dir_out, paste0(spp_info$srvy, "_Test_Species/"))
+  if (dir.exists(test_dir)) unlink(test_dir, recursive = TRUE)
+
+  expect_no_error({
+    resample_tests(
+      spp_dfs = spp_dfs,
+      spp_info = spp_info,
+      grid_yrs = grid_yrs,
+      dir_out = dir_out,
+      test = TRUE,
+      parallel = FALSE,
+      n_knots = NULL,
+      model_type = "wrapper_sdmtmb"
+    )
+  })
+
+  # Check output files exist
+  expect_true(file.exists(file.path(test_dir, "fit_df.csv")))
+  expect_true(file.exists(file.path(test_dir, "fit_pars.csv")))
+  expect_true(file.exists(file.path(test_dir, "fit_check.csv")))
+  expect_true(file.exists(file.path(test_dir, "index.csv")))
+
+  # Check that parquet files exist
+  expect_true(file.exists(file.path(test_dir, "df_1.parquet")))
+  expect_true(file.exists(file.path(test_dir, "df_2.parquet")))
+
+  # Read and check CSV contents
+  fit_df <- read.csv(file.path(test_dir, "fit_df.csv"))
+  expect_true("estimate" %in% names(fit_df))
+  fit_index <- read.csv(file.path(test_dir, "index.csv"))
+  expect_true(all(c("year", "index") %in% names(fit_index)))
+})
+
+test_that("select_n_knots returns expected values", {
+  expect_equal(select_n_knots(500), 50)
+  expect_equal(select_n_knots(3000), 100)
+  expect_equal(select_n_knots(9000), 200)
+  expect_equal(select_n_knots(20000), 500)
+  expect_equal(select_n_knots(100000), 1000)
+})
+
diff --git a/vignettes/a-simple-example.Rmd b/vignettes/a-simple-example.Rmd
@@ -245,7 +245,16 @@ The `purrr::map` function is used to apply the `clean_and_resample` function to
 The number you input for `n_knots` can make or break your model. We have more details on considerations when choosing `n_knots` or allowing the function to select the number of knots for you in the [Importance of `n_knots` in `{sdmTMB}` Models section](#importance-of-sdmtmb-models). 
 :::
 
-```{r run-models, eval = FALSE}
+You can run this in testing mode
+```{r}
+clean_and_resample(spp_info = spp_list,
+                   catch,
+                   seq_from, seq_to, seq_by, 
+                   tot_dataframes, replicate_num,
+                   grid_yrs, dir_out, test = TRUE, n_knots = 300) 
+```
+
+```{r run-models}
 start.time <- Sys.time()
 purrr::map(
   seq_len(nrow(spp_list)),
@@ -256,19 +265,11 @@ purrr::map(
     n_knots = 300
   )
 )
-write.csv(
-  x = data.frame(
-    time = as.numeric(Sys.time() - start.time),
-    units = units(Sys.time() - start.time)
-  ),
-  file = paste0(dir_final, srvy, "_simple_time.csv")
-)
 ```
 
 Don't freak out if this takes a long time. It should take the following amount of time:
-```{r time}
-a <- read.csv(file = paste0(dir_final, srvy, "_simple_time.csv"))
-print(paste0("Completed in: ", round(a$time, 2), " ", a$units))
+```{r, echo = FALSE}
+print(paste0("Completed in: ", round(as.numeric(Sys.time() - start.time), 2), " ", units(Sys.time() - start.time)))
 ```
 
 ```{r sink-results-backup, include = FALSE}
diff --git a/vignettes/outputEBS_walleye_pollock/fit_check.csv b/vignettes/outputEBS_walleye_pollock/fit_check.csv
@@ -0,0 +1,4 @@
+"srvy","common_name","model_fn","model_family","model_anisotropy","model_spatiotemporal","effort","hessian_ok","eigen_values_ok","nlminb_ok","range_ok","gradients_ok","se_magnitude_ok","se_na_ok","sigmas_ok","all_ok"
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_1",FALSE,TRUE,TRUE,TRUE,TRUE,FALSE,TRUE,TRUE,FALSE
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_2",FALSE,TRUE,TRUE,TRUE,TRUE,FALSE,TRUE,TRUE,FALSE
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_3",FALSE,TRUE,TRUE,TRUE,TRUE,FALSE,TRUE,TRUE,FALSE
diff --git a/vignettes/outputEBS_walleye_pollock/fit_df.csv b/vignettes/outputEBS_walleye_pollock/fit_df.csv
@@ -0,0 +1,13 @@
+"srvy","common_name","model_fn","model_family","model_anisotropy","model_spatiotemporal","effort","term","estimate","std.error","conf.low","conf.high"
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_1","factor(year)2021",20.3035148497297,1020.3623946318,-1979.57002980765,2020.17705950711
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_1","factor(year)2022",21.4673087246129,2860.24471327961,-5584.50931627451,5627.44393372374
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_1","factor(year)2023",21.006938931454,2166.32132286837,-4224.90483283172,4266.91871069463
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_1","factor(year)2024",20.7459415141075,1707.14795575973,-3325.20256805615,3366.69445108436
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_2","factor(year)2021",21.7210323464771,2442.61938463255,-4765.72498947272,4809.16705416567
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_2","factor(year)2022",21.6242224589619,NA,NA,NA
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_2","factor(year)2023",21.5654413409925,NA,NA,NA
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_2","factor(year)2024",21.4913033210022,1669.97594763063,-3251.60140908318,3294.58401572518
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_3","factor(year)2021",20.4135309519096,1679.53172014375,-3271.40815142244,3312.23521332626
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_3","factor(year)2022",20.4193816640038,1170.4479263295,-2273.61639972141,2314.45516304942
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_3","factor(year)2023",20.4135309512272,NA,NA,NA
+"EBS","walleye pollock","total_catch_wt_kg ~ 0 + factor(year)","delta_gamma",TRUE,"iid, iid","05_3","factor(year)2024",20.3390234786106,1783.02587428651,-3474.327473626,3515.00552058322
diff --git a/vignettes/outputEBS_walleye_pollock/fit_pars.csv b/vignettes/outputEBS_walleye_pollock/fit_pars.csv
diff --git a/vignettes/outputEBS_walleye_pollock/index.csv b/vignettes/outputEBS_walleye_pollock/index.csv
