FLR

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

library(FLCandy)

FishLife

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

Plotting and Data wrangling

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

library(ggplot2)
theme_set(theme_bw())

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

library(fishblicc)

Back to Top

par=lhPar(FLPar(linf=100))
eq=lhEql(par)

om=as(eq,"FLStock")
om=window(fwd(om,fbar=m(om)[4,-1],sr=eq),start=19)

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

dimnames(om)$year =an(dimnames(om)$year)+1980
range(om)["plusgroup"]=range(om)["maxage"]
om=fwd(om,fbar=window(fbar(om),start=2021)%*%
                 FLQuant(seq(1,10,length.out=61)),sr=eq)
dim(window(om,end=2080))

[1] 41 82  1  1  1  1

ak =invALK(par,cv=0.1,age=an(dimnames(om)$age),bin=1)
lfd=lenSamples(catch.n(om),ak,n=5000)

units(lfd)="cm"

#plotLengths(group(lfd, sum, year=year-year%%10))

dat=subset(transform(as.data.frame(lfd),lustrum=5*(year%/%5)),data>0) 

dt2=ddply(dat,.(lustrum), with, {
  lmodeFn<-function(len,n,bins=25) {
    
    dat=data.frame(bin=cut(len,breaks=bins),n=n)
    res=ddply(dat,.(bin), with, data.frame(freq=sum(n)))
    res=as.character(subset(res,freq==max(freq))[1,"bin"])
    unbin(res)$mid}

  unbin=function(x){
  left =as.numeric(substr(x,2,unlist(gregexpr(",",x))-1))
  right=as.numeric(substr(x,  unlist(gregexpr(",",x))+1,nchar(x)-1))
  mid  =(left+right)/2
  
  data.frame(left=left,right=right,mid=mid,n=x)}

  lc=lmodeFn(len,data)*0.5
  data.frame(lc   =lc,
             lmean=weighted.mean(len,data*(len>=lc)))})

gghistogram(dat,x="len",weight="data",bins=120)+
  coord_flip()+
  scale_x_reverse()+
  facet_grid(.~lustrum,scale="free")+xlab("Length (cm)")+
  geom_vline(aes(xintercept=lc),    data=dt2,col="red")+
  geom_vline(aes(xintercept=lmean), data=dt2,col="blue")+
  #scale_x_continuous(limits=c(0,60))+
  theme_bw(16)+
  theme(legend.position = "none", 
                    axis.title.x = element_blank(), 
                    axis.text.x  = element_blank(), 
                    axis.ticks.x = element_blank())

library(remotes)
install_github("PaulAHMedley/fishblicc")

library(fishblicc)

## Prepare some data in the required format
frq=subset(lfd,year==2000&iter==1&len%in%10:101)$data

frq=c(34,42,45,68,65,56,91,123,126,126,113,123,107,111,113,103,91,85,74,53,43,
      48,35,38,26,31,21,25,20,17,22,16,19,18,19,12,11,11,20,11,7,11,7,13,15,6,9,
      7,8,4,7,6,8,3,3,5,6,7,5,3,7,7,2,5,1,3,2,2,0,0,1,2,0,1,3,2,3,2,2,2,0,1,2,2,
      1,0,1,2,0,0,0)
  
dl=blicc_dat(
  #model_name ="Base: Mk ~ Length-inverse",
  LLB        = 10:100,           # Lower boundaries of length bins
  fq         = list(`Test`=frq), # frequency data
  Linf       = c(100,10),
  sel_fun    = c("logistic"),    # Selectivity functions
  Catch      = c(1)              # Relative catch for each gear
  #gear_names = c("test"),
  #Mk         = 1.5,# Natural mortality prior mean (at reference length)
  #ref_length =50,  # Reference length for the length-inverse natural mortality model
  #a          =0.3, # Length-weight scale parameter (optional)
  #b          =3,   # Length-weight exponent
  #L50        =52   # Length at 50% maturity
)

## Fit the model to these data 
slim <- blicc_mpd(dl)

rp_res <- blicc_ref_pts(slim, dl)

plot_expected_frequency(rp_res, gear=1)


stf=blicc_fit(dl, ntarget=100, nwarmup=200, nchain=1)

Back to Top

Author information

Laurence Kell. laurie@seaplusplus.co.uk

Acknowledgements

Software Versions

R version 4.2.1 (2022-06-23)

• knitr: 1.43

• FLCore: 2.6.19.9105

• FLBRP: 2.5.9.9022

• FLasher: 0.7.1.9221

• FLife: 3.4.0

• ggplotFL: 2.7.0.9132

• FLCandy: 0.1.0

• rfishbase: 3.1.9.99

• SPMpriors: 1.0.4

• FishLife: 3.0.0

• plyr: 1.8.9

• dplyr: 1.1.4

• reshape: 0.8.9

• ggplot2: 3.4.4

• ggpubr: 0.6.0

• ggcorrplot: 0.1.3

• GGally: 2.1.2

Compiled: Sun Jan 21 12:52:01 2024

Back to Top