FAO: Small Pelagics

01

Jump to Introduction

Jump to Packages

Jump to Operating Model Conditioning

Jump to Indicators

Jump to Classification Skill

Jump to More Information

Jump to References

Introduction

This vignette, conditions an Operating Model based on life history relationships, where the parameters are extracted from FishBase using the FishLife package. Although, the example is for Rastrineobola argentea the code is generic will work for any species. However, the process requires user inout before the Operating Model should be used to provide advice.

First an Operating Model, based on an age structured FLStock object is created. Then an Observation Error Model simulates pseudo data, for catch (\(C_t\)), an index of relative abundance (\(I_t\)), a length based-indicator (LBI) and relative harvest rate (i.e. \(H_t = C_t/I_t\)).

The simulated time series are compared to the “true” values from the Operating Model, to identify which indicators can track changes in status, i.e. by plotting the

• relative index of abundance verses vulnerable biomass (\(V_t\))
• length based-indicator \(L_{mean}\) v Fishing Mortality (\(F_t\))
• relative harvest rate verses \(F_t\)

Back to Top

Packages

setwd("~/Desktop/active/fao/Rmd/")

FLR packages

library(FLCore)
library(FLBRP)
library(FLasher)
library(FLife)
library(ggplotFL)

library(FLCandy)

Packages for extracting life history parameters from FishLife

library(rfishbase)
library(SPMpriors)
library(FishLife)

Utility packages, for data wrangling

library(plyr)
library(dplyr)
library(reshape)

Plotting packages

library(ggplot2)
theme_set(theme_bw())

library(ggpubr)
library(ggcorrplot)
library(GGally)

Installation

Libraries (i.e. packages) can be installed from cran or via gihtub using the remotes package for the FLR, FishBase and FishLife packages

install.packages("remotes")

library(remotes)

remotes::install_github("flr/FLCore")

There are compatability problems with FishBase, so install an earlier version

remotes::install_github("ropensci/rfishbase",ref="d150f2e0f5")

remotes::install_github("james-thorson/FishLife")

Operating Model Conditioning

Life History Parameters

Read in life history parameters from Fishbase using the FishLife package

pars=flmvn_traits(Genus="Rastrineobola",Species="argentea",Plot=FALSE)

     [,1]           [,2]         
[1,] "K"            "M"          
[2,] "Winfinity"    "Loo"        
[3,] "tmax"         "tm"         
[4,] "Lm"           "Temperature"
[5,] "ln_margsd"    "rho"        
[6,] "logitbound_h" "ln_r"       

Extract the life history parameters required to parameterise the Operating Model, these are parameters for the individual growth (\(k\) and \(L_{\infty}\)), age at maturity (a_{50}\(, natural mortality (\)M\(), the steepness of the stock recruitment relationship (\)h$), and recruitment variability (sigR).

par=FLPar(c(cast(pars[[1]][,c("trait","mu.sp")],
                 ~trait)[,c("K","Loo","tm","M","h","sigR")]))

Rename parameters, so they are compatible with the FLR package naming convention, and then save the expected values, and the covariances

dimnames(par)[[1]]=c("k","linf","a50","m","s","sigmaR")
cov=pars[[2]][[1]][[1]][c("K","Loo","tm","M","h"),c("K","Loo","tm","M","h")]
dimnames(cov)=list(c("k","linf","a50","m","s"),
                   c("k","linf","a50","m","s"))
R.argentea=list(param=par,cov=cov)

par =R.argentea[[1]]
par =rbind(lhPar(rbind(par[c("k","linf","a50","s")],FLPar("t0"=-0.1))),par["m"])

par["a50"]=1

par

An object of class "FLPar"
params
     linf         k        t0         a         b       l50       a50     ato95 
  14.4809    0.7467   -0.1000    0.0003    3.0000    6.9898    1.0000    1.0000 
     asym        bg        m1        m2        m3         s         v      sel1 
   1.0000    3.0000    0.5500   -1.6100    1.4400    0.9843 1000.0000    1.7827 
     sel2      sel3         m 
   1.0000 5000.0000    1.5870 
units:  NA 

Natural mortality rates (\(M\)) vary with body size and age, therefore, take the single value from FishBase and use the ‘generalized length-inverse mortality’ (GLIM, Lorenzen 2022) paradigm to model M-at-age, i.e.

\(log(M)= m_1 + m_2 \times log(L/L_{50}_{\infty}) + m_3 \times log(k)\)

par["m2"]=-1.30
par["m3"]= 1.08

Scale m1 so that \(M\) at \(l_{50}\) is equal to the ‘FishBase’ value

$ m_1 = log(M) - m_2log(L_{50}/L_{}) - m_3log(k)$

par["m1"]=log(par["m"])%-%(par["m2"]%*%log(par["l50"]%/%par["linf"]))%-%(par["m3"]%*%log(par["k"]))

par

An object of class "FLPar"
params
     linf         k        t0         a         b       l50       a50     ato95 
  14.4809    0.7467   -0.1000    0.0003    3.0000    6.9898    1.0000    1.0000 
     asym        bg        m1        m2        m3         s         v      sel1 
   1.0000    3.0000   -0.1696   -1.3000    1.0800    0.9843 1000.0000    1.7827 
     sel2      sel3         m 
   1.0000 5000.0000    1.5870 
units:  NA 

ages  =FLQuant(0:10,dimnames=list(age=0:10))
length=vonB(ages,par)

ggplot(as.data.frame(FLQuants(length=vonB(ages,par),
                              weight=len2wt(length,par),
                              mat   =logistic(length,par),
                              m     =mGlim(length,par))))+
  geom_line(aes(age,data))+
  facet_wrap(~qname,scale="free",ncol=2)

Figure 1 Vectors of quantities-at-age.

ggplot(as.data.frame(FLQuants(mat        =logistic(length,par),
                              Selectivity=dnormal(ages,par))))+
  geom_line(aes(age,data))+
  facet_wrap(~qname,scale="free",ncol=2)

Figure 2 Selectivity and maturity-at-age.

Monte Carlo Simulations

nits=101

par=R.argentea[[1]]
cov=R.argentea[[2]]

save(par,cov,nits,file="~/Desktop/tmp/t.RData")
load("~/Desktop/tmp/t.RData")

par=rlmvn(nits,par[dimnames(cov)[[1]]],diag(cov)^0.5*c(0.1,0,0,0,0),cov2cor(cov))

par=rbind(lhPar(rbind(par[c("k","linf","a50","s")],propagate(FLPar("t0"=-0.1),nits))),par["m"])

par["m1"]= 0.28
par["m2"]=-1.30
par["m3"]= 1.08
par["m1"]=log(par["m"])%-%(par["m2"]%*%log(par["l50"]%/%par["linf"]))%-%(par["m3"]%*%log(par["k"]))

eq=lhEql(par,m=mGlimFn) 

plot(eq,refpts="msy",labels=FALSE,panels=3)

Figure 3 Production functions.

fbar(eq)[]=par["m"]
om=as(brp(eq),"FLStock")
om=fwd(om,f=fbar(om)[,-1],sr=eq)

trnd=FLQuant(c(t(maply(c(seq(2.5,1.0,length.out=51),
                         seq(1,  0.5,length.out=51)[-1]), function(x)
  seq(1,x,length.out=21)))),dimnames=list(year=seq(21),iter=seq(101)))

kick=FLQuant(rep(c(trnd[,21]),each=21),dimnames=list(year=seq(21),iter=seq(101)))

oms=FLStocks("Constant"=om)
oms[["Trend"]]   =fwd(om,f=fbar(om)[,ac(51:71)]%*%trnd,sr=eq)
oms[["Kick"]]    =fwd(om,f=fbar(om)[,ac(51:71)]%*%kick,sr=eq)
oms=FLStocks(llply(oms, window, start=31, end=71))

[1] "maxfbar has been changed to accomodate new plusgroup"
[1] "maxfbar has been changed to accomodate new plusgroup"
[1] "maxfbar has been changed to accomodate new plusgroup"

p1=plot(window(oms[-1],start=41),metrics=list( 
        #LBI  =function(x) 1/fbar(x),
         F    =fbar,
         Stock=ssb),iter=c(25,51,75))+
    theme_bw(16)+
    theme(legend.position="none")+
    scale_fill_manual(  "Scenario",values=rainbow(3)[-2])+
    scale_colour_manual("Scenario",values=c(rainbow(3)[-2],c("grey50","black","grey50")))+
  facet_grid(qname~stock,scale="free_y")
p1  

Figure 4 Operating Model, for the changes in fishing mortality.

Back to Top

Indicators

source("~/Desktop/flr/FLCandy/R/saveTmp.R")
sv(oms)
ld()


x=oms[[2]]
Index=ssb(x)%*%rlnorm(dim(x)[6],iter(vb(x),1)%=%0,0.3)

p2=plot(window(oms[-1],start=41),metrics=list(
          LBI =function(x) {res=indicators.len(x[-1], indicators=c('lmean'), 
                                          params=apply(par,1,median), metric='stock.n')[[1]] 
                            names(dimnames(res))[1]="age"
                            res},
          Index=function(x) rlnorm(dim(x)[6],log(ssb(x)),0.3))
         #,HR   =function(x) apply(catch(x),c(2,6),sum)%/%ssb(x))
         ,iter=c(25,51,75),alpha=c(0.1,0.1))+
    theme_bw(16)+ 
    theme(legend.position="none")+
    scale_fill_manual(  "Scenario",values=rainbow(3)[-2])+
    scale_colour_manual("Scenario",values=c(rainbow(3)[-2],c("grey50","black","grey50")))+
    scale_y_continuous(limit=c(9.0,NA)) 
p2+facet_grid(qname~stock,scale="free_y") 

Figure 5 Indices generated by the Observation Error Model, for the three trends in fishing mortality.

Back to Top

Classification Skill

pars=apply(par,1,median)

ts=ldply(oms[c("Kick","Trend")], function(x){
   
  index=rlnorm(dim(x)[6],log(vb(x)),0.3)
  
  rbind.fill(
     cbind("What"="LBI",model.frame(FLQuants(x,
          OM =function(x) fbar(x),
          OEM=function(x) {res=indicators.len(x[-1], indicators=c('lmean'), 
                            params=pars, metric=catch.n)[[1]] 
                      names(dimnames(res))[1]="age"
                      1/res}),drop=T)),
     cbind("What"="Index",model.frame(FLQuants(
          OM =vb(x),
          OEM=index),drop=T)),
     cbind("What"="HR",model.frame(FLQuants(x,
          OM =function(x) fbar(x),
          OEM=function(x) catch(x)%/%index),drop=T)))})


ts =ddply(ts,.(.id,What), transform, OM =(OM -mean(OM)) /var(OM)^0.5, 
                                     OEM=(OEM-mean(OEM))/var(OEM)^0.5)
ggplot(aes(OM,OEM),data=subset(ts,What!="HR"))+  
  geom_point()+
  facet_wrap(What~.id,scale="free",ncol=2)+
  geom_smooth(col="blue",method=lm,se=FALSE)+
  geom_hline(aes(yintercept=0),col="red",linetype=2)+
  geom_vline(aes(xintercept=0),col="red",linetype=2)+
  theme_bw(16)

sv(ts)
ld() 
roc=ddply(ts,.(.id,What,year), with, roc(OM+1,OEM+1))

ggplot(subset(roc,What!="HR"&year%in%51:55))+  
  geom_path(aes(FPR,TPR,col=ac(year),group=ac(year)))+
  facet_grid(What~.id)+
  geom_abline(aes(slope=1,intercept=0),col="red")+
  theme_bw(24)+
  theme(legend.position="bottom")+
  scale_colour_manual("Year",values=rainbow(6))

Figure 6 ROC curves

Back to Top

More Information

• You can submit bug reports, questions or suggestions on FLife at the FLife issue page [^1], or on the FLR mailing list.
• Or send a pull request to https://github.com/lauriekell/FLife/
• For more information on the FLR Project for Quantitative Fisheries Science in R, visit the FLR webpage [^2].
• The latest version of FLife can always be installed using the devtools package, by calling

Author information

Laurence Kell. laurie@seaplusplus.co.uk

Acknowledgements

Software Versions

• R version 4.2.1 (2022-06-23) R: r version$version.string

knitr: r packageVersion(‘knitr’)

FLCore: r packageVersion(‘FLCore’) FLBRP: r packageVersion(‘FLBRP’) FLasher: r packageVersion(‘FLasher’) FLife: r packageVersion(‘FLife’) ggplotFL: r packageVersion(‘ggplotFL’) FLCandy: r packageVersion(‘FLCandy’)

rfishbase: r packageVersion(‘rfishbase’) SPMpriors: r packageVersion(‘SPMpriors’) FishLife: r packageVersion(‘FishLife’)

plyr: r packageVersion(‘plyr’) dplyr: r packageVersion(‘dplyr’) reshape: r packageVersion(‘reshape’)

ggplot2: r packageVersion(‘ggplot2’) ggpubr: r packageVersion(‘ggpubr’) ggcorrplot: r packageVersion(‘ggcorrplot’) GGally: r packageVersion(‘GGally’)

Compiled: r format(Sys.time(), ‘%d %B, %Y’)

• Compiled: Fri Jan 12 12:39:51 2024

Back to Top

References

Back to Top

Lorenzen, Kai. 2022. “Size-and Age-Dependent Natural Mortality in Fish Populations: Biology, Models, Implications, and a Generalized Length-Inverse Mortality Paradigm.” Fisheries Research 255: 106454.