Running a retrospective analysis on a FIMS model
retrospective.RmdSet up FIMS model
Following the minimal demo vignette, we will set up a FIMS model and run it.
# Load sample data
data("data1")
# Prepare data for FIMS model
data_4_model <- FIMSFrame(data1)
# Create parameters
parameters <- data_4_model |>
create_default_configurations() |>
create_default_parameters(data = data_4_model)
# Run the model with optimization
base_model <- parameters |>
initialize_fims(data = data_4_model) |>
fit_fims(optimize = TRUE)
#> ✔ Starting optimization ...
#> ℹ Restarting optimizer 3 times to improve gradient.
#> ℹ Maximum gradient went from 0.02469 to 0.00037 after 3 steps.
#> ✔ Finished optimization
#> ✔ Finished sdreport
#> ℹ FIMS model version: 0.7.1
#> ℹ Total run time was 9.99432 seconds
#> ℹ Number of parameters: fixed_effects=77, random_effects=0, and total=77
#> ℹ Maximum gradient= 0.00037
#> ℹ Negative log likelihood (NLL):
#> • Marginal NLL= 3240.1261
#> • Total NLL= 3240.1261
#> ℹ Terminal SB= 1803.27712
# Clear memory post-run
clear()
retro_fit <- run_fims_retrospective(years_to_remove = 0:2,
data_4_model, parameters,
n_cores = 3)
#> ...Running in parallel with multicore
#> ℹ running model with 0 years of data removed
#> ✔ Starting optimization ...
#> ℹ Restarting optimizer 3 times to improve gradient.
#> ℹ Maximum gradient went from 0.02469 to 0.00037 after 3 steps.
#> ✔ Finished optimization
#> ✔ Finished sdreport
#> ℹ FIMS model version: 0.7.1
#> ℹ Total run time was 10.28581 seconds
#> ℹ Number of parameters: fixed_effects=77, random_effects=0, and total=77
#> ℹ Maximum gradient= 0.00037
#> ℹ Negative log likelihood (NLL):
#> • Marginal NLL= 3240.1261
#> • Total NLL= 3240.1261
#> ℹ Terminal SB= 1803.27712
#> ℹ running model with 1 years of data removed
#>
#> ✔ Starting optimization ...
#> ℹ Restarting optimizer 3 times to improve gradient.
#> ℹ Maximum gradient went from 0.03303 to 0.00161 after 3 steps.
#> ✔ Finished optimization
#> ✔ Finished sdreport
#> ℹ FIMS model version: 0.7.1
#> ℹ Total run time was 26.60289 seconds
#> ℹ Number of parameters: fixed_effects=77, random_effects=0, and total=77
#> ℹ Maximum gradient= 0.00161
#> ℹ Negative log likelihood (NLL):
#> • Marginal NLL= 3140.71461
#> • Total NLL= 3140.71461
#> ℹ Terminal SB= 2315.0298
#> ℹ running model with 2 years of data removed
#>
#> ✔ Starting optimization ...
#> ℹ Restarting optimizer 3 times to improve gradient.
#> ℹ Maximum gradient went from 0.01823 to 0.00303 after 3 steps.
#> ✔ Finished optimization
#> ✔ Finished sdreport
#> ℹ FIMS model version: 0.7.1
#> ℹ Total run time was 27.97789 seconds
#> ℹ Number of parameters: fixed_effects=77, random_effects=0, and total=77
#> ℹ Maximum gradient= 0.00303
#> ℹ Negative log likelihood (NLL):
#> • Marginal NLL= 3034.31988
#> • Total NLL= 3034.31988
#> ℹ Terminal SB= 3115.03058
# plot the SSB time series from each model
# NOTE: this doesn't work because we don't have a year value.
# there's a time column but it's empty
library(ggplot2)
retro_fit[["estimates"]] |>
dplyr::filter(label == "spawning_biomass") |>
dplyr::select(label, year_i, estimated, retro_year) |>
dplyr::group_by(label, retro_year) |>
dplyr::filter(dplyr::row_number() <= dplyr::n() - retro_year) |>
dplyr::ungroup() |>
ggplot(aes(x = year_i, y = estimated, group = as.factor(retro_year))) +
geom_line(aes(color = as.factor(retro_year)), linewidth = 1.2) +
stockplotr::theme_noaa(discrete = T)
# get SSB time series from each model
# proof that values are different among models
# retro_fit[["estimates"]] |>
# dplyr::filter(label == "spawning_biomass") |>
# dplyr::select(label, year_i, estimated, retro_year) |>
# tidyr::pivot_wider(names_from = retro_year, values_from = estimated) |>
# tail()
# SSBtable
# label year_i `0` `1` `2`
# <chr> <int> <dbl> <dbl> <dbl>
# 1 spawning_biomass 26 2001. 1960. 2213.
# 2 spawning_biomass 27 1939. 1889. 2201.
# 3 spawning_biomass 28 1914. 1854. 2232.
# 4 spawning_biomass 29 1798. 1725. 2153.
# 5 spawning_biomass 30 1901. 1815. 2623.
# 6 spawning_biomass 31 1803. 2315. 3115.