gadget-framework / rgadget

Rgadget is a set of useful utilities for gadget
6 stars 12 forks source link

Gadget projections #99

Open bthe opened 4 years ago

bthe commented 4 years ago

I've started the redevelopment of the forward projections functionality with the more flexible gadget_project_* functions. The logic behind these functions is that the user defines the projections, unlike gadget.forward, in a stepwise manner thus giving greater control over individual processes:

All these functions could be applied repeatedly if the number of stocks in the model is greater than 2. @vbartolino and @pamelajwoods would you mind giving these functions a go and report on any shortcomings.

vbartolino commented 4 years ago

when recruitment is scheduled in a timestep different from step 1 the following error is returned:

Error: Can't subset columns that don't exist.
✖ The columns `venimm.rec.pre.2019.1`, `venimm.rec.pre.2020.1`, `venimm.rec.pre.2021.1`, `venimm.rec.pre.2022.1`, `venimm.rec.pre.2023.1`, etc. don't exist.
Run `rlang::last_error()` to see where the error occurred.

I think the issue is that gadget_project_recruitment still looks for projecting recruitment in the first step rather than the correct one

dplyr::mutate(year = sprintf('%s.rec.pre.%s.1',stock,.data$year))
bthe commented 4 years ago

Actually the recruitment step is not changed, meaning that if the immature stock was recruited at step 3 it will continue to recruit at step 3. These functions, however, only assume one recruitment step.

The error you have looks like you haven't yet called gadget on the new model variant, so directly after gadget_project_fleets you need to call gadget with gadget_evaluate to get the update parameter file:

gadget_project_time() %>% 
  gadget_project_stocks(imm.file = 'gssimm',mat.file = 'gssmat') %>% 
  gadget_project_fleets() %>% 
  gadget_evaluate(params.out = paste(attr(.,'variant_dir'),'params.pre',sep='/'),
                  params.in = 'WGTS/params.final')

ps: stock venimm sounds pretty ominous.

vbartolino commented 4 years ago

this is the code I'm using (borrowed from the example) and that returns the error

res <- 
      gadget_project_time() %>% 
      gadget_project_stocks(imm.file='venimm', mat.file='venmat') %>% 
      gadget_project_fleets(pre_fleet='comven2', fleet_type="linearfleet") %>% 
      gadget_evaluate(params.out = paste(attr(.,'variant_dir'),'params.pre',sep='/'),
                      params.in = 'WGTS/params.final') %>% 
      gadget_project_recruitment(stock = 'venimm', 
                                 recruitment = fit$stock.recruitment,
                                 params.file = paste(attr(.,'variant_dir'),'params.pre',sep='/'))

If I check the file params.pre which is created along the way I can see that the recruitment parameters

...
venimm.rec.pre.2019.2              1       -9999     9999        0
venimm.rec.pre.2020.2              1       -9999     9999        0
venimm.rec.pre.2021.2              1       -9999     9999        0
...

recruitment happens in step 2 in this model like correctly reported by these parameters but judging from the error the function is looking for venimm.rec.pre.2019.1`, `venimm.rec.pre.2020.1, .... If I hack the function gadget_projection_recruitment and change the line dplyr::mutate(year = sprintf('%s.rec.pre.%s.1',stock,.data$year)) with dplyr::mutate(year = sprintf('%s.rec.pre.%s.2',stock,.data$year)) at least this error disappear. Although I notice that the format of the file params.pre is completely changed...

bthe commented 4 years ago

The most recent version I've fixed this issue in gadget_project_recruitment among a few other brainfarts.

vbartolino commented 4 years ago

The fix does the job, great! Few more to report. Based on the code below (adapted from your example):

res <- 
      gadget_project_time(num_years=15) %>% 
      gadget_project_stocks(imm.file='venimm', mat.file='venmat') %>% 
      gadget_project_fleets(pre_fleet='comven2', fleet_type="linearfleet") %>% 
      gadget_evaluate(params.out = paste(attr(.,'variant_dir'),'params.pre',sep='/'),
                      params.in = 'WGTS/params.final') %>% 
      gadget_project_recruitment(stock = 'venimm', n_replicates=100,
                                 recruitment = fit$stock.recruitment,
                                 params.file = paste(attr(.,'variant_dir'),'params.pre',sep='/')) %>%
      gadget_project_advice(harvest_rate = 1:10/100, n_replicates=100,
                       pre_fleet = "comven2",
                       params.file = paste(attr(.,'variant_dir'),'params.pre',sep='/')) %>%
      gadget_project_output(imm.file='venimm',mat.file='venmat', pre_fleets="comven2") %>% 
      gadget_evaluate(params.in = paste(attr(.,'variant_dir'),'params.pre',sep='/'))  %>%
      read.printfiles(paste(attr(.,'variant_dir'),'out',sep='/'))

I need to specify the pre_fleets argument in gadget_project_output otherwise the gadget_evaluate that follows doesn't work and return the following message

Warning message:
In system(run.string, intern = TRUE, ignore.stderr = ignore.stderr) :
  running command '/usr/local/R350/lib/R/library/gadget2/bin/gadget  -s    -i prj/params.pre  -main prj/main          2>/dev/null' had status 1

The code above run until the second last line. It prints the files under 'out' and at a quick check they look fine, but when I try to read them with read.printfiles I get the following error:

Error in !suppress : invalid argument type
In addition: Warning message:
In .f(.x[[i]], ...) :
  Old style printfile detected, update Gadget to the most recent version

which is surprising given that I've just updated gadget using the installation via the R package gadget2 - very handy by the way ;)

bthe commented 4 years ago

Yes, this is my bad. Note that the notation x %>% f(y=g(.)) is equivalent to x %>% f(a=., y=g(.)). If you don't want x to be considered as the first argument of the function call you need to wrap the call in curly braces i.e like this x %>% {f(y=g(.))}. In our case we need to wrap the call to read.printfiles in curly braces:

... %>%
{read.printfiles(paste(attr(.,'variant_dir'),'out',sep='/'))} %>% 
... 

it will otherwise think that the output from the preceding steps is the output directory and the output directory is supplied into the suppress argument. Hope this fixes the issue.

vbartolino commented 4 years ago

yes, it fix it, thanks! Interesting use of {...}, I didn't think about it. I need to do some more tests, but so far so good and it really improves forecasts functionalities.

bthe commented 4 years ago

I've added a functionality to gadget_project_fleets that allows the user to define a common multiplier to the fleet data:

... %>% 
gadget_project_fleets(pre_fleet = 'fleet1',
                        common_mult = 'effort.pre', ## defining a common multiplier for (all) fleets
                        pre_propotion = 0.5, ## fleet specific scalar
                        ...) %>%
... %>%

You can then just use the gadget_project_advice like the group of fleets is a single fleet:

... %>%
  gadget_project_advice(pre_fleet = 'effort',
                        harvest_rate = 0.2,
                        params.file = paste(attr(.,'variant_dir'),'params.pre',sep='/'),
                        advice_cv = 0.2) %>%
...