This vignette includes an example of how to specify metapopulation dynamics in multi-wham. Now the stocks are able to move to another region and contribute to the SSB during the spawning season (i.e. reproductive mixing).
library(wham)
library(whamMSE)
main.dir = here::here()
info <- generate_basic_info(base.years = 1:50, move_dyn = 1) # now stocks don't need to move back during the spawning season
basic_info = info$basic_info # collect basic information
catch_info = info$catch_info # collect fleet catch information
index_info = info$index_info # collect survey information
F_info = info$F # collect fishing information
# see more details using ?generate_basic_info
basic_info <- generate_NAA_where(basic_info = basic_info, move.type = 2) # "bidirectional" movement
move <- generate_move(basic_info = basic_info, move.type = 2, move.rate = c(0.2,0.1), move.re = "constant")
n_stocks <- as.integer(basic_info['n_stocks'])
n_regions <- as.integer(basic_info['n_regions'])
n_fleets <- as.integer(basic_info['n_fleets'])
n_indices <- as.integer(basic_info['n_indices'])
n_ages <- as.integer(basic_info['n_ages'])
sigma <- "rec+1"
re_cor <- "iid"
ini.opt <- "equilibrium"
sigma_vals <- array(0.05, dim = c(n_stocks, n_regions, n_ages)) # Assume a low NAA sigma
# Set para. for B-H function
alpha <- 12
beta <- 1.5e-4
NAA_re <- list(N1_model=rep(ini.opt,n_stocks),
sigma=rep(sigma,n_stocks),
cor=rep(re_cor,n_stocks),
recruit_model = 3,
recruit_pars = rep(list(c(alpha,beta)),n_stocks), # assume same B-H s-r functions for all stocks
sigma_vals = sigma_vals) # NAA_where must be specified in basic_info!
input <- prepare_wham_input(basic_info = basic_info, NAA_re = NAA_re, move = move,
catch_info = catch_info, index_info = index_info, F = F_info)
random = input$random # check what processes are random effects
input$random = NULL # so inner optimization won't change simulated RE
om <- fit_wham(input, do.fit = F, do.brps = T, MakeADFun.silent = TRUE)
# Note: do.fit must be FALSE (no modeling fitting yet)
set.seed(123)
input <- om$input
input$data <- om$simulate(complete=TRUE)
plot(input$data$SSB[,1], type = "l", col = "red",
ylab = "SSB", xlab = "Years",
main = "SSB", ylim = c(0,max(input$data$SSB)))
lines(input$data$SSB[,2], type = "l", col = "blue")
plot(input$data$NAA[1,1,,1], type = "l", col = "red",
ylab = "Recruitment", xlab = "Years",
main = "Recruitment", ylim = c(0,max(input$data$NAA)))
lines(input$data$NAA[2,2,,1], type = "l", col = "blue")