lkerr / groundfish-MSE

Development of Robust Management Strategies for Northeast Groundfish Fisheries in a Changing Climate
MIT License
4 stars 2 forks source link

CAA non-functioning #201

Open mle2718 opened 2 years ago

mle2718 commented 2 years ago

I can't get a CAA to run. Not sure who to assign, so I assigned it to alot of people.

Part 1 on troubleshooting:

  1. Modified set_om_parameters_global.R to read in mprocTest.csv
    mprocfile<-"mprocTest.csv"
  2. Modified runSim to pick up lines 3 and 4 (the CAA models)
    mproc_bak<-mproc
    mproc<-mproc_bak[3:4,]
  3. Ran runSim.R from inside Rstudio. Here is my error message:
    Error in if ((y - fmyearIdx)%%mproc[m, "AssessFreq"] == 0) { : 
    argument is of length zero

    traceback gives:

    get_advice(stock = stock[[i]])

    This makes sense: the management procedures table needs an AssessFreq column.

Part 2:

I was hoping that creating an AssessFreq column in the management procedures file would fix the problem. I was not successful.

  1. Restart R.
  2. Created a modified version of mprocTest.csv called mprocTest2.csv. Added a column with AssessFreq = 2 for all rows.
  3. Modified ``set_om_parameters_global.R'' to read in mprocTest2.csv
  4. Ran runSim.R from inside Rstudio. Here is my error message:
    Error in .Call("getParameterOrder", data, parameters, new.env(), NULL,  : 
    Incorrect number of arguments (4), expecting 3 for 'getParameterOrder'
    Called from: MakeADFun(data = data, parameters = parameters, DLL = "caa", 
    map = map_par, silent = TRUE)

    and traceback gives:

    7: MakeADFun(data = data, parameters = parameters, DLL = "caa", 
       map = map_par, silent = TRUE) at get_caa.R#51
    6: eval(substitute(expr), e)
    5: eval(substitute(expr), e)
    4: within.list(stock, {
       invisible(capture.output(dyn.load(dynlib("assessment/caa"))))
       map_par <- list(log_M = factor(NA), log_oe_sumCW = factor(NA), 
           log_oe_sumIN = factor(NA), log_pe_R = factor(NA))
       lb <- tmb_lb
       ub <- tmb_ub
       active_idx <- !names(lb) %in% names(map_par)
       for (i in which(active_idx)) {
           tmb_par[[i]] <- get_svNoise(tmb_par[[i]], cv = startCV, 
               lb = tmb_lb[[i]], ub = tmb_ub[[i]])
       }
       data <- tmb_dat
       parameters <- tmb_par
       obj <- MakeADFun(data = data, parameters = parameters, DLL = "caa", 
           map = map_par, silent = TRUE)
       opt <- try(nlminb(start = obj$par, objective = obj$fn, gradient = obj$gr, 
           lower = unlist(lb[active_idx]), upper = unlist(ub[active_idx]), 
           control = list(iter.max = 1e+05, eval.max = 1000)))
       sdrep <- try(sdreport(obj))
       if (class(sdrep)[1] == "try-error") 
           browser()
       rep <- obj$report()
       rep$gradient.fixed <- sdrep$gradient.fixed
       rep$maxGrad <- max(abs(rep$gradient.fixed))
       rep_data <- rep[c("obs_sumCW", "obs_paaCN", "obs_sumIN", 
           "obs_paaIN", "obs_effort", "laa", "waa", "slxI", "timeI")]
       rep_par <- rep[c("log_M", "R_dev", "log_ipop_mean", "ipop_dev", 
           "oe_sumCW", "oe_sumIN", "log_pe_R", "log_selC", "log_qC", 
           "log_qI")]
       rep_NLL <- rep[c("NLL", "NLL_sumCW", "NLL_sumIN", "NLL_R_dev", 
           "NLL_paaCN", "NLL_paaIN")]
       d2 <- data.frame(rep$sumCW, rep$sumIN, rep$paaCN, rep$paaIN)
       if (platform != "Linux") {
           Sys.setenv(PATH = path0)
       }
    })
    3: within(stock, {
       invisible(capture.output(dyn.load(dynlib("assessment/caa"))))
       map_par <- list(log_M = factor(NA), log_oe_sumCW = factor(NA), 
           log_oe_sumIN = factor(NA), log_pe_R = factor(NA))
       lb <- tmb_lb
       ub <- tmb_ub
       active_idx <- !names(lb) %in% names(map_par)
       for (i in which(active_idx)) {
           tmb_par[[i]] <- get_svNoise(tmb_par[[i]], cv = startCV, 
               lb = tmb_lb[[i]], ub = tmb_ub[[i]])
       }
       data <- tmb_dat
       parameters <- tmb_par
       obj <- MakeADFun(data = data, parameters = parameters, DLL = "caa", 
           map = map_par, silent = TRUE)
       opt <- try(nlminb(start = obj$par, objective = obj$fn, gradient = obj$gr, 
           lower = unlist(lb[active_idx]), upper = unlist(ub[active_idx]), 
           control = list(iter.max = 1e+05, eval.max = 1000)))
       sdrep <- try(sdreport(obj))
       if (class(sdrep)[1] == "try-error") 
           browser()
       rep <- obj$report()
       rep$gradient.fixed <- sdrep$gradient.fixed
       rep$maxGrad <- max(abs(rep$gradient.fixed))
       rep_data <- rep[c("obs_sumCW", "obs_paaCN", "obs_sumIN", 
           "obs_paaIN", "obs_effort", "laa", "waa", "slxI", "timeI")]
       rep_par <- rep[c("log_M", "R_dev", "log_ipop_mean", "ipop_dev", 
           "oe_sumCW", "oe_sumIN", "log_pe_R", "log_selC", "log_qC", 
           "log_qI")]
       rep_NLL <- rep[c("NLL", "NLL_sumCW", "NLL_sumIN", "NLL_R_dev", 
           "NLL_paaCN", "NLL_paaIN")]
       d2 <- data.frame(rep$sumCW, rep$sumIN, rep$paaCN, rep$paaIN)
       if (platform != "Linux") {
           Sys.setenv(PATH = path0)
       }
    }) at get_caa.R#3
    2: get_caa(stock = tempStock) at get_advice.R#10
    1: get_advice(stock = stock[[i]])

Part 3:

I noticed there were more columns in the mproc.csv.I added those to mprocTest3.csv, make the appropriate change to ``set_om_parameters_global.R'' and tried to run runSim.R again. Here are the error messages:

Error in .Call("getParameterOrder", data, parameters, new.env(), NULL,  : 
Incorrect number of arguments (4), expecting 3 for 'getParameterOrder'

Called from: MakeADFun(data = data, parameters = parameters, DLL = "caa", 
    map = map_par, silent = TRUE)

> traceback()
7: MakeADFun(data = data, parameters = parameters, DLL = "caa", 
       map = map_par, silent = TRUE) at get_caa.R#51
6: eval(substitute(expr), e)
5: eval(substitute(expr), e)
4: within.list(stock, {
       invisible(capture.output(dyn.load(dynlib("assessment/caa"))))
       map_par <- list(log_M = factor(NA), log_oe_sumCW = factor(NA), 
           log_oe_sumIN = factor(NA), log_pe_R = factor(NA))
       lb <- tmb_lb
       ub <- tmb_ub
       active_idx <- !names(lb) %in% names(map_par)
       for (i in which(active_idx)) {
           tmb_par[[i]] <- get_svNoise(tmb_par[[i]], cv = startCV, 
               lb = tmb_lb[[i]], ub = tmb_ub[[i]])
       }
       data <- tmb_dat
       parameters <- tmb_par
       obj <- MakeADFun(data = data, parameters = parameters, DLL = "caa", 
           map = map_par, silent = TRUE)
       opt <- try(nlminb(start = obj$par, objective = obj$fn, gradient = obj$gr, 
           lower = unlist(lb[active_idx]), upper = unlist(ub[active_idx]), 
           control = list(iter.max = 1e+05, eval.max = 1000)))
       sdrep <- try(sdreport(obj))
       if (class(sdrep)[1] == "try-error") 
           browser()
       rep <- obj$report()
       rep$gradient.fixed <- sdrep$gradient.fixed
       rep$maxGrad <- max(abs(rep$gradient.fixed))
       rep_data <- rep[c("obs_sumCW", "obs_paaCN", "obs_sumIN", 
           "obs_paaIN", "obs_effort", "laa", "waa", "slxI", "timeI")]
       rep_par <- rep[c("log_M", "R_dev", "log_ipop_mean", "ipop_dev", 
           "oe_sumCW", "oe_sumIN", "log_pe_R", "log_selC", "log_qC", 
           "log_qI")]
       rep_NLL <- rep[c("NLL", "NLL_sumCW", "NLL_sumIN", "NLL_R_dev", 
           "NLL_paaCN", "NLL_paaIN")]
       d2 <- data.frame(rep$sumCW, rep$sumIN, rep$paaCN, rep$paaIN)
       if (platform != "Linux") {
           Sys.setenv(PATH = path0)
       }
   })
3: within(stock, {
       invisible(capture.output(dyn.load(dynlib("assessment/caa"))))
       map_par <- list(log_M = factor(NA), log_oe_sumCW = factor(NA), 
           log_oe_sumIN = factor(NA), log_pe_R = factor(NA))
       lb <- tmb_lb
       ub <- tmb_ub
       active_idx <- !names(lb) %in% names(map_par)
       for (i in which(active_idx)) {
           tmb_par[[i]] <- get_svNoise(tmb_par[[i]], cv = startCV, 
               lb = tmb_lb[[i]], ub = tmb_ub[[i]])
       }
       data <- tmb_dat
       parameters <- tmb_par
       obj <- MakeADFun(data = data, parameters = parameters, DLL = "caa", 
           map = map_par, silent = TRUE)
       opt <- try(nlminb(start = obj$par, objective = obj$fn, gradient = obj$gr, 
           lower = unlist(lb[active_idx]), upper = unlist(ub[active_idx]), 
           control = list(iter.max = 1e+05, eval.max = 1000)))
       sdrep <- try(sdreport(obj))
       if (class(sdrep)[1] == "try-error") 
           browser()
       rep <- obj$report()
       rep$gradient.fixed <- sdrep$gradient.fixed
       rep$maxGrad <- max(abs(rep$gradient.fixed))
       rep_data <- rep[c("obs_sumCW", "obs_paaCN", "obs_sumIN", 
           "obs_paaIN", "obs_effort", "laa", "waa", "slxI", "timeI")]
       rep_par <- rep[c("log_M", "R_dev", "log_ipop_mean", "ipop_dev", 
           "oe_sumCW", "oe_sumIN", "log_pe_R", "log_selC", "log_qC", 
           "log_qI")]
       rep_NLL <- rep[c("NLL", "NLL_sumCW", "NLL_sumIN", "NLL_R_dev", 
           "NLL_paaCN", "NLL_paaIN")]
       d2 <- data.frame(rep$sumCW, rep$sumIN, rep$paaCN, rep$paaIN)
       if (platform != "Linux") {
           Sys.setenv(PATH = path0)
       }
   }) at get_caa.R#3
2: get_caa(stock = tempStock) at get_advice.R#10
1: get_advice(stock = stock[[i]])

Here's the output of sessionInfo()

R version 4.0.5 (2021-03-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Red Hat Enterprise Linux 8.5 (Ootpa)

Matrix products: default
BLAS/LAPACK: /usr/lib64/libopenblasp-r0.3.12.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8    LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] fishmethods_1.11-3  fGarch_3042.83.2    fBasics_3042.89.1   timeSeries_3062.100 timeDate_3043.102   ASAPplots_0.2.18   
 [7] data.table_1.14.2   forcats_0.5.1       stringr_1.4.0       dplyr_1.0.7         purrr_0.3.4         readr_2.1.1        
[13] tidyr_1.1.4         tibble_3.1.6        ggplot2_3.3.5       tidyverse_1.3.1     glue_1.6.0          abind_1.4-5        
[19] TMB_1.8.0           tmvtnorm_1.5        gmm_1.6-6           sandwich_3.0-1      Matrix_1.3-4        mvtnorm_1.1-3      
[25] msm_1.6.9          

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.7          lubridate_1.8.0     lattice_0.20-45     zoo_1.8-9           assertthat_0.2.1    utf8_1.2.2         
 [7] R6_2.5.1            cellranger_1.1.0    backports_1.4.0     reprex_2.0.1        httr_1.4.2          pillar_1.6.4       
[13] rlang_0.4.12        spatial_7.3-14      readxl_1.3.1        minqa_1.2.4         rstudioapi_0.13     nloptr_1.2.2.3     
[19] splines_4.0.5       lme4_1.1-27.1       munsell_0.5.0       broom_0.7.10        numDeriv_2016.8-1.1 compiler_4.0.5     
[25] modelr_0.1.8        pkgconfig_2.0.3     tidyselect_1.1.1    expm_0.999-6        fansi_0.5.0         crayon_1.4.2       
[31] tzdb_0.2.0          dbplyr_2.1.1        withr_2.4.3         MASS_7.3-54         grid_4.0.5          nlme_3.1-153       
[37] jsonlite_1.7.2      gtable_0.3.0        lifecycle_1.0.1     DBI_1.1.1           magrittr_2.0.1      scales_1.1.1       
[43] cli_3.1.0           stringi_1.7.6       fs_1.5.1            xml2_1.3.3          ellipsis_0.3.2      generics_0.1.1     
[49] vctrs_0.3.8         boot_1.3-28         tools_4.0.5         hms_1.1.1           survival_3.2-13     colorspace_2.0-2   
[55] bootstrap_2019.6    rvest_1.0.2         haven_2.4.3 

You can look it over in https://github.com/lkerr/groundfish-MSE/tree/CAA_non_functional

lkerr commented 2 years ago

Some googling of this issue revealed that it may have to do with the matrix and TBM package versions (e.g. if Matrix version is more recent than the most recently released TMB version).

Have you gotten a single run to work on your local machine (not linux server)? Wondering if you get this same error in a local run?

samtruesdell commented 2 years ago

Hi Min-Yang -- any luck on this? @mle2718

mle2718 commented 2 years ago

Not yet, although I just gave it a first shot this am.

  1. I don't think it's the Matrix and TMB problem that @lkerr described. I have TMB 1.8.0, which is from 2022-03-07 and Matrix 1.3.-4, which is from 2021-05-24. So TMB is more recent than Matrix.

  2. I noticed that my TMB 1.8.0 was built in R 4.1.2 even though the server is running R4.0.5. I have no idea how this happened. But I also tried installing TMB 1.8.1, both to check if there was something fixed in 1.8.1 or if the R version mismatches were causing the problem. No luck there.

  3. I tried switching servers, but that's a non-starter because it's an older R version.

  4. I'm also using an older R version on my laptop -- I'll try to get a newer R version installed and see if that fixes the problem.

If that doesn't work, I really don't know what to do: I'm waay out of my element here in troubleshooting this.

samtruesdell commented 2 years ago

OK if a newer R on your laptop doesn't work I'll give it a spin and see if I can figure out what's going on

mle2718 commented 2 years ago

Identical error message on my windows laptop. Here's the sessionInfo()

R version 4.2.0 (2022-04-22 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19044)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.utf8  LC_CTYPE=English_United States.utf8    LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C                           LC_TIME=English_United States.utf8    

attached base packages:
[1] stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] fishmethods_1.11-3  fGarch_3042.83.2    fBasics_3042.89.1   timeSeries_3062.100 timeDate_3043.102   ASAPplots_0.2.18   
 [7] data.table_1.14.2   forcats_0.5.1       stringr_1.4.0       dplyr_1.0.8         purrr_0.3.4         readr_2.1.2        
[13] tidyr_1.2.0         tibble_3.1.6        ggplot2_3.3.5       tidyverse_1.3.1     glue_1.6.2          abind_1.4-5        
[19] TMB_1.8.1           tmvtnorm_1.5        gmm_1.6-6           sandwich_3.0-1      Matrix_1.4-1        mvtnorm_1.1-3      
[25] msm_1.6.9          

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.8.3        lubridate_1.8.0     lattice_0.20-45     zoo_1.8-10          assertthat_0.2.1    utf8_1.2.2         
 [7] R6_2.5.1            cellranger_1.1.0    backports_1.4.1     reprex_2.0.1        httr_1.4.2          pillar_1.7.0       
[13] rlang_1.0.2         spatial_7.3-15      readxl_1.4.0        minqa_1.2.4         rstudioapi_0.13     nloptr_2.0.0       
[19] splines_4.2.0       lme4_1.1-29         munsell_0.5.0       broom_0.8.0         numDeriv_2016.8-1.1 compiler_4.2.0     
[25] modelr_0.1.8        pkgconfig_2.0.3     tidyselect_1.1.2    expm_0.999-6        fansi_1.0.3         crayon_1.5.1       
[31] tzdb_0.3.0          dbplyr_2.1.1        withr_2.5.0         MASS_7.3-56         grid_4.2.0          nlme_3.1-157       
[37] jsonlite_1.8.0      gtable_0.3.0        lifecycle_1.0.1     DBI_1.1.2           magrittr_2.0.3      scales_1.2.0       
[43] cli_3.2.0           stringi_1.7.6       fs_1.5.2            xml2_1.3.3          ellipsis_0.3.2      generics_0.1.2     
[49] vctrs_0.4.1         boot_1.3-28         tools_4.2.0         hms_1.1.1           survival_3.3-1      colorspace_2.0-3   
[55] bootstrap_2019.6    rvest_1.0.2         haven_2.5.0      
samtruesdell commented 2 years ago

I will try to give this a look. Are you working on a branch that is pushed?

mle2718 commented 2 years ago

CAA_non_functional. But that is just master with 1-2 very small changes to make it a working example.

samtruesdell commented 2 years ago

Hi Min-Yang -- I did some work on this a few days ago and was stymied; I cannot get the TMB-based CAA model to work and haven't been able to track down the problem. @mdmazur or @ahart1 -- did either of you use this version of the CAA model for simulations? @mle2718 is it possible for you to use the ASAP option? Slower, but I'm sure that's what's been used lately.

mle2718 commented 2 years ago

Ugh. I feel a little better that I'm not the only one with the problem. But not much better.

This is ASAP, right?

samtruesdell commented 2 years ago

Yes that's the ASAP. If you choose assessment model class "ASAP" in mproc rather than "CAA" it will fit that model. "CAA" is pretty much the same structurally, just generally simpler so it runs (ran!) faster. I was mostly using CAA, but the last running branch that uses it might be quite a bit behind. Mackenzie or Amanda can jump in here, but I think WINE was required for the way things were set up, but that might just be because the available compilation was windows and I would imagine there are other compiled versions available.

lkerr commented 2 years ago

ASAP should work for you @mle2718 , but will likely be too slow unless you are running on a cluster. Let me know if this is not working for you with the run time and we will do a bit more investigating on CAA and circle back.

lkerr commented 2 years ago

One more point to do with versioning of packages, we have a compilation of the Rlibraries used when running on the HPCC (ie. consistent versions of all packages that Mackenzie was using in her runs). Let me know if you want us to transfer these to you.

mle2718 commented 2 years ago

Tim Miller has set me up with the source code for ASAP3. He also said this:

If the group wants to try using WHAM let me know. WHAM can take an ASAP dat file and fit a comparable model with maybe 2 lines of R code.

samtruesdell commented 2 years ago

Awesome -- WHAM might be worthwhile if there are any gains in speed. I think you can just substitute within functions/assessment/get_ASAP to pull from the Linux compiled executable instead of using wine.

mle2718 commented 2 years ago

From TIm:

wham (using TMB) isn't any faster than ADMB for models without random effects. TMB has a big speed up over ADMB when mixed effects models are estimated.

mdmazur commented 2 years ago

Sorry about the delay on this, I was out at sea. I haven't been able to find out what the problem is either. I get an error when sourcing runSim but not when in browser() mode. The CAA model was working when I originally merged things into the master but I've made some changes since then that could be part of the issue. The old version of the master is the branch CopyofOldMaster. It does not have all the new features, but the CAA model should work. I can also try to find the push that broke the CAA, but that might take sometime.

mle2718 commented 2 years ago

Thanks @mdmazur. Now that ASAP is running, a non-functioning CAA model is not a big deal for me. Of course, other people are going to have different priorities.

mle2718 commented 1 year ago

I'll see if I can force install an "older" matrix version. I can also try running it locally to troubleshoot.

But in the long run, it's not a great solution. The economics module just takes too long to run locally.

On Friday, April 15, 2022, Kerr Fisheries Research Lab < @.***> wrote:

Some googling of this issue revealed that it may have to do with the matrix and TBM package versions (e.g. if Matrix version is more recent than the most recently released TMB version).

Have you gotten a single run to work on your local machine (not linux server)? Wondering if you get this same error in a local run?

— Reply to this email directly, view it on GitHub https://github.com/lkerr/groundfish-MSE/issues/201#issuecomment-1100068271, or unsubscribe https://github.com/notifications/unsubscribe-auth/AIT5KA447KYW6TCHE7T4UJDVFFMDZANCNFSM5TLKBZXQ . You are receiving this because you were assigned.Message ID: @.***>

-- Sent from Gmail Mobile

samtruesdell commented 1 year ago

Does this indicate you are planning to go back to CAA? There are some vagaries there that make it a bit different from ASAP (and more similar to the Great Lakes models I was working with before).