Skip to content

Problem with catch scenarios for longer management periods #154

@alko989

Description

@alko989

There are sometimes issues when the management period is longer than one year and Catch fractiles are used. The F/Fmsy at the end of the management period is higher than 1.

It seems to affect assessments with low biomass in the end of the time series.

A possible solution is to split Cp and dtcp to yearly chunks instead of a single long catch duration.

ReprEx

library(spict)
#> Loading required package: TMB
#> Welcome to spict_v1.3.5@8a4e1b
set.seed(1234)
nt <- 50
inp <- list(nseasons=4, splineorder=3)
s1 <- seq(0, 10, by=1)
s2 <- seq(11, nt-3-1/inp$nseasons, by=1/inp$nseasons)
s3 <- seq(nt-3, nt-1, by = 1)
inp$timeC <- c(s1,s2,s3)
inp$dtc <- c(rep(1, length(s1)),
             rep(0.25, length(s2)),
             rep(1, length(s3)))
inp$timeI <- seq(0.1, nt-1/inp$nseasons, by=0.5)
inp$ini <- list(logK=log(1000), logm=log(123), logq=log(1), logn=log(2),
                logbkfrac=log(0.7), logsdf=log(0.123), logF0=log(0.4321),
                logphi=log(c(0.3, 0.5, 1.8)))
inp$dteuler <- 1/8
inpsim <- sim.spict(inp)
#> The stabilise option was used for fitting / is specified in the input list (inp$stabilise). This activates six very vague priors. Acknolowdging these priors when simulating is not yet implemented. 
#> Additional priors were used for fitting / are specified in the input list (inp$priors): logn, logalpha, logbeta. Acknolowdging these priors when simulating is not yet implemented.
inpsim$dtc <- inp$dtc ## Very important, sim.spict drops it
inpsim$maneval <- nt+5
inpsim$maninterval <- c(nt,nt+5)

inpsim$priors$logbkfrac <- c(log(0.7), 0.001, 1)
inpsim$priors$logn <- c(log(2), 0.001, 1)
inpsim$dtc <- inp$dtc
inpsim <- check.inp(inpsim)

plotspict.data(inpsim)
#> Error in xy.coords(x, y, xlabel, ylabel, log): 'x' and 'y' lengths differ

fit <- fit.spict(inpsim)
fit <- add.man.scenario(fit, "Fmsy")
fit <- add.man.scenario(fit, "C35", fractiles = list(catch = 0.35))
fit <- add.man.scenario(fit, "All35", fractiles = list(catch = 0.35, ffmsy=0.35, bbmsy=0.35))
fit <- add.man.scenario(fit, "FB35_noC", fractiles = list( ffmsy=0.35, bbmsy=0.35))
plot2(fit) ## This does not, indxmax is not defined

sumspict.manage(fit)
#> SPiCT timeline:
#>                                                   
#>       Observations              Management        
#>      0.00 - 50.00            50.00 - 55.00        
#>  |-----------------------| ----------------------|
#> 
#> Management evaluation: 55.00
#> 
#> Predicted catch for management period and states at management evaluation time:
#> 
#>                     C B/Bmsy F/Fmsy
#> 1. Fish at Fmsy 166.1   0.27   1.00
#> 2. C35          148.1   0.02   2.53
#> 3. All35        145.3   0.02   2.66
#> 4. FB35_noC     163.3   0.30   0.92

Created on 2022-03-08 by the reprex package (v2.0.1)

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions