library(knitr)
## Global options
opts_chunk$set(cache     =TRUE,
               cache.path='cache/turbot-om/',
               echo      =TRUE,
               eval      =TRUE,
               prompt    =FALSE,
               comment   =NA,
               message   =FALSE,
               warning   =FALSE,
               tidy      =FALSE,
               fig.height=6,
               fig.width =8,
               fig.path  ='tex/turbot-om-')

iFig=0

Equilibrium Dynamics

Time Series

Stochasticity

References

Introduction

In Management Strategy Evaluation (MSE) a simulation model, known as an Operating Model (OM), is used to represent the resource dynamics and the uncertainty about them. The OM is then used to simulate pseudo data to allow Harvest Control Rules (HCRs) and Management Procedures (MPs) to be simulation tested using closed loop simulation. A MP is the combination of monitoring data, analysis method and management measures as well as reference points and harvest control rules. Linking the OM and the MP requires an Observation Error Model (OEM) to generate fishery-dependent or fishery-independent resource monitoring data. The OEM reflects the uncertainties, between the actual dynamics of the resource and perceptions arising from observations and assumptions by modelling the differences between for example of an index of abundance and the actual value in the OM.

The OM can also be used for open loop simulation testing, where psuedo data are simulated and used by an assessment procedure to estimate stock status but not fed back into the OM.

FLife

Here an OM is conditioned on turbot life history parameters using the FLife package.

library(FLife)
turbot=FLPar(c(linf=59.1,  k=0.28, t0=-0.4, a=0.01111,
               b   =3.15,a50=4.0, l50=43.25),units="NA")

The lhPar method is then used to derive the parameters for natural mortality-at=age, based on @gislason2008does, and default parameters and relationships for selection pattern and stock recruitment.

turbot=lhPar(turbot)

turbot

The default parameters can be changed, e.g. by changing a parameter. sl is the standard deviation for the lefthand limb of the double normal selection pattern, here we change it from 2 to 1 to make it steeper.

Back to Top

Equilibrium Dynamics {#Equilibrium}

http://127.0.0.1:13112/rmd_output/4/#Time The parameters are then used by lhEql to simulate the equilibrium dynamics by combining the spawner/yield per recruit relationships with a stock recruiment relationship.

library(FLBRP)
eq=lhEql(turbot)

range(eq)[c("minfbar","maxfbar")]=ceiling(mean(turbot["sel1"]))
library(ggplotFL)
sel<-function(x) 
  catch.sel(x)%/%fapex(catch.sel(x))

ggplot(FLQuants(eq,"catch.wt","mat","m","catch.sel"=sel))+
  geom_line(aes(age,data))+
  facet_wrap(~qname,scale="free")+
  scale_x_continuous(limits=c(0,20))+ 
  guides(colour=guide_legend(title="Species",title.position="top"))

Figure r iFig=iFig+1; iFig Vectors of m, selection pattern, maturity and weight-at-age.

Estimate equilibrium dynamics and reference points

plot(eq,refpts="msy")

Figure r iFig=iFig+1; iFig Expected, equilibrium, dynamics and reference points.

Back to Top

Time Series {#Time}

To go from equilibrium to time series dynamics the FLBRP object created by lhEql can be coerced to an FLStock object.

First change the F time series so that it represents a time series where the stock was origionally lightly exploited, F increased until the stock was overfished and then fishing pressure was reduced to ensure spawning stock biomass was greater than $B_{MSY}$.

library(FLasher)
fbar(eq)=refpts(eq)["msy","harvest"]%*%FLQuant(c(rep(.1,19),
                                              seq(.1,2,length.out = 30),
                                              seq(2,1.0,length.out = 10),
                                              rep(1.0,61)))[,1:105]
om=as(eq,"FLStock")

om=fwd(om,fbar=fbar(om)[,-1],sr=eq)

Back to Top

Stochasticity {#Stochasticity}

To simulation random variation in the time series, deviations around the stock recruitment relationship was modelled as a random variable.

These deviates were then used to create a stochastic time series by projecting the dynamics from year 1.

set.seed(1234)
srDev=rlnoise(100,fbar(om)%=%0,.3,b=0.0)
plot(srDev,iter=c(7,2,9))

Figure r iFig=iFig+1; iFig Time series of recruitment deviates

om=propagate(om,dim(srDev)[6])
om=fwd(om,fbar=fbar(om)[,-1],deviances=srDev,sr=eq)
plot(om,iter=1:3)

Figure r iFig=iFig+1; iFig Stochastic Time series of F, SSB, recruitment and yield

While to generate data for use in the MP, random measurement error was added to the simulated catch per unit effort (CPUE).

save(om,file="../data/om.RData",compress="xz")
save(eq,file="../data/eq.RData",compress="xz")

Back to Top

\newpage

Software Versions

Author information

Laurence Kell. laurie@seaplusplus.es

Acknowledgements

This vignette and many of the methods documented in it were developed under the MyDas project funded by the Irish exchequer and EMFF 2014-2020. The overall aim of MyDas is to develop and test a range of assessment models and methods to establish Maximum Sustainable Yield (MSY) reference points (or proxy MSY reference points) across the spectrum of data-limited stocks.

References {#References}

Back to Top



flr/mydas documentation built on Jan. 19, 2024, 10:33 a.m.