LBSPR: Length-Based Spawning Potential Ratio

Package Description

This package is well described in http://barefootecologist.com.au/lbspr.html

The authors of LBSPR state: “The LBSPR method has been developed for data-limited fisheries, where few data are available other than a representative sample of the size structure of the vulnerable portion of the population (i.e., the catch) and an understanding of the life history of the species. The LBSPR method does not require knowledge of the natural mortality rate (M), but instead uses the ratio of natural mortality and the von Bertalanffy growth coefficient (K) (M/K), which is believed to vary less across stocks and species than M (Prince et al. 2015).”

Spawning Potential Ratio (SPR) is the ratio of total reproductive production at equilibrium for a given level of fishing mortality divided by the reproductive production in the unfished state.

The authors of LBSPR note that “The stock-recruitment relationship and the steepness parameter are difficult to estimate and unlikely to be known with any certainty for data-poor stocks. Several studies have explored the levels of SPR to be used as reference points, and it is generally accepted that a SPR of 35-40% is sustainable for most species (Clark, 2002, 1993; Mace and Sissenwine, 1993), although a risk-adverse target SPR may be higher for species thought to have low steepness (e.g., Dorn, 2002). For our example, the SPR target (SPRtarg) was set to 0.40, which corresponds to a %SSB of ~30-40% for the most commonly observed values of steepness (h = 0.6; Clark, 2002; Myers et al., 2002).”

This package is run from R as found on CRAN. R can be installed https://cran.r-project.org/.

The recommended way to run R is using RStudio. After you have installed R, then install RStudio. It can be found https://www.rstudio.com/products/rstudio/download/.

Installation

There are two ways to run LBSPR - a) either on the Barefoot Ecologist web site as a shiny app http://barefootecologist.com.au/lbspr (scroll to the bottom of the page) or b) within R which means you need to install the LBSPR package. These instructions are for running the test data set within R.

NB Once you have installed LBSPR, you don’t need to run the install.packages code chunk again.

Note: Version for R4.0.0 and above is not yet available for this package !

Check_Install.packages <- function(pkg){
  new.pkg <- pkg[!(pkg %in% installed.packages()[, "Package"])]
  if (length(new.pkg)) 
    install.packages(new.pkg,dependencies = TRUE) else print(paste0("'",pkg,"' has been installed already!"))
  sapply(pkg, require, character.only = TRUE)
}

Check_Install.packages("LBSPR")
## Installing package into 'C:/~/R/win-library/4.0'
## (as 'lib' is unspecified)
## package 'LBSPR' successfully unpacked and MD5 sums checked
## 
## The downloaded binary packages are in
##  C:\Users\den123\AppData\Local\Temp\Rtmpia1htn\downloaded_packages
## Loading required package: LBSPR
## Warning: package 'LBSPR' was built under R version 4.0.3
## LBSPR 
##  TRUE
library(LBSPR)

On the CRAN website (which describes the R package version on https://cran.r-project.org/web/packages/LBSPR/vignettes/LBSPR.html), the package can be used to a) simulate the expected length composition, growth curve, SPR and yield curves using the LBSPR model and b) fit to empirical length data to provide an estimate of the spawning potential ratio (SPR) to simulate expected life-history length composition.

The description below is to fit our empirical length data using LBSPR. Although we describe the steps below, there is a good description in https://cran.r-project.org/web/packages/LBSPR/vignettes/LBSPR.html#2_first_steps, which we recommend if you want to run LBSPR’s own test data that is downloaded with the package.

Reading in the data

The following are required to read in the data:
1. Create the blank LB_pars object

MyPars <- new("LB_pars")
## A blank LB_pars object created
## Default values have been set for some parameters
slotNames(MyPars) # shows the elements or slots in LB_pars
##  [1] "Species"      "MK"           "M"            "Linf"         "L_units"     
##  [6] "CVLinf"       "L50"          "L95"          "Walpha"       "Walpha_units"
## [11] "Wbeta"        "FecB"         "Steepness"    "Mpow"         "R0"          
## [16] "SL50"         "SL95"         "MLL"          "sdLegal"      "fDisc"       
## [21] "FM"           "SPR"          "BinMin"       "BinMax"       "BinWidth"

Although MYPars has 25 slots, which are all named, only 6 of these are needed to fit the model.

  1. Fill in the life history parameters.

Our dataset was created using Stock Synthesis based on the following life-history parameters. You can put your own values here, but those below are ours. LBSPR is a female-only model so we are input the female values for the biologIcal parameter. We also need to set the Target Reference Point is usually about 30-40% as described above. We used 40% in this example.

MyPars@Species <- "OurCodoid" # Our own special species ;-)
MyPars@Linf <- 81.44505 # Asymptotic length
MyPars@L50 <- 45 # Length at which 50% of the females are mature
MyPars@L95 <- 56.77776 # Length at which 95% of the females are mature
MyPars@MK <- 1.333333333 # Natural mortality/von Bertalnffy growth <em>K</em> coefficial ratio
MyPars@SPR <- 0.4 # The target reference point
MyPars@L_units <- "cm"
  1. Creating a blank length object and read in your length data

You need to work out if you have:
+ a single year or multiple years of length-frequency data, + length-frequencies or raw lengths, and + whether the data has a header or not.

For length-frequency data, the first column of the data set must contain the mid-points of the length bins, and the remaining columns contain the counts for each length class.

If you have raw data, then you must specify the minimum and maximum in sizes and the bin widths which are the final three slots in MyPars. Each of these uses different code which can be found in Section 4.2: https://cran.r-project.org/web/packages/LBSPR/vignettes/LBSPR.html#42_reading_in_example_csv

The data type (frequencies or raw) must be specified in the call to the new function.

A valid LB_pars object must also be supplied to the new function.

Following the normal convention in R, if a header row is present in the data file, the argument header in the call to new must be set to TRUE.

The format of the LBSPR examples are described here:

In our case, we have multiple (31) years of length-frequency data (with headers). This dataset needs first to have its first four columns (“year”,“Fleet”,“Sex”,“Stage1_wght”) removed, and then to be transposed to have the first column contain the size bins and then the next 31 columns contain 31 years of length-frequency data.

NB We are assuming that you have saved our test data in the same directory as the R markdown code and that you have opened a new R session. Otherwise you can change our code to point to the directory where the data file is and call it your working directory using the setwd command.

temp <- read.csv("LengthData.csv",header=T) # Reading in the length-frequency data
LenFreq <- as.matrix(temp[,5:21])
mids <- seq(12.5, by = 5, length.out = ncol(LenFreq))
LenFreq <- cbind(mids,t(LenFreq))
colnames(LenFreq)<- c(NA,1971:2001)
colnames(LenFreq)<-as.numeric(colnames(LenFreq))

# If your length frequency data is in the correct format, you can use the same file naming protocol as given in online examples.
Len2 <- new("LB_lengths", LB_pars=MyPars, LenFreq, dataType="freq", header=TRUE)
## Warning in if (class(file) == "character") {: the condition has length > 1 and
## only the first element will be used

Fitting the model to our data

Here are some considerations when fitting the model.

The default model type (modtype) is Growth-Type_Group (GTG), which uses Growth-Type_Groups to account for size-based selectivity (Hordyik et al. 2015c). This model also has the ability to include variable natural mortality at size. This model typically estimates a lower fishing mortality for a given size-structure compared to the age-structured model.

The other option is a conventional age-structured equilibrium population model (absel). An important assumption of this model structure is that selectivity is age- rather than length-based (Hordyk et al. 2015a, b).

To see all the model options, the help function below will send you to the help web links:

help("LBSPRfit")
## starting httpd help server ... done
help("LBSPRfit_")

We now fit the model:

plotSize(Len2) # Plots the length data if you want
## Warning: Years must be numeric values
## Attempting to convert to numeric values

myFit2 <- LBSPRfit(MyPars, Len2) # Fits the model, here we use all the defaults
## Warning: Years must be numeric values
## Attempting to convert to numeric values
## Fitting model
## Year:
## 1
## 2
## 3
## 4
## 5
## 6
## 7
## 8
## 9
## 10
## 11
## 12
## 13
## 14
## 15
## 16
## 17
## 18
## 19
## 20
## 21
## 22
## 23
## 24
## 25
## 26
## 27
## 28
## 29
## 30
## 31

Summarise and plot the results

myFit2@Ests # smoother parameter estimates
##        SL50  SL95   FM  SPR
##  [1,] 45.14 61.86 0.26 0.78
##  [2,] 45.21 62.24 0.27 0.78
##  [3,] 44.72 61.46 0.28 0.78
##  [4,] 44.05 59.89 0.29 0.77
##  [5,] 44.18 59.49 0.34 0.74
##  [6,] 44.82 60.03 0.41 0.70
##  [7,] 45.70 61.22 0.48 0.67
##  [8,] 46.50 62.50 0.54 0.64
##  [9,] 46.93 63.27 0.58 0.63
## [10,] 47.00 63.54 0.61 0.62
## [11,] 46.99 63.69 0.65 0.61
## [12,] 46.67 63.42 0.71 0.59
## [13,] 45.92 62.26 0.78 0.57
## [14,] 45.37 61.23 0.90 0.54
## [15,] 45.19 60.88 1.07 0.49
## [16,] 45.14 60.95 1.28 0.45
## [17,] 44.86 60.82 1.51 0.41
## [18,] 44.05 59.77 1.75 0.37
## [19,] 43.10 58.24 2.05 0.33
## [20,] 42.39 56.16 2.41 0.30
## [21,] 42.94 56.10 2.91 0.27
## [22,] 44.09 57.03 3.45 0.25
## [23,] 45.29 58.26 3.95 0.25
## [24,] 46.19 59.36 4.40 0.25
## [25,] 44.95 57.05 3.98 0.25
## [26,] 45.44 57.27 3.91 0.24
## [27,] 47.00 59.13 4.03 0.24
## [28,] 49.02 61.81 4.17 0.25
## [29,] 51.09 64.77 4.29 0.26
## [30,] 52.79 67.33 4.32 0.28
## [31,] 54.03 69.31 4.38 0.29
data.frame(rawSL50=myFit2@SL50, rawSL95=myFit2@SL95, rawFM=myFit2@FM, rawSPR=myFit2@SPR) # individual point estimates for each year
##    rawSL50 rawSL95 rawFM    rawSPR
## 1    44.87   58.64  0.15 0.8394017
## 2    50.77   73.89  0.31 0.7513949
## 3    46.61   69.38  0.16 0.8352910
## 4    36.07   48.20  0.00 0.9999976
## 5    39.02   50.05  0.17 0.7898145
## 6    42.39   53.46  0.36 0.6534303
## 7    46.55   60.46  0.60 0.5680123
## 8    50.20   67.54  0.73 0.5601673
## 9    50.47   68.25  0.67 0.5828741
## 10   47.88   64.67  0.42 0.6662015
## 11   49.98   68.00  0.55 0.6249850
## 12   51.04   72.26  0.59 0.6162004
## 13   43.98   60.99  0.29 0.7120676
## 14   41.61   54.48  0.38 0.6326437
## 15   43.86   56.64  0.70 0.4983426
## 16   47.53   62.93  1.09 0.4285089
## 17   50.16   70.04  1.37 0.4020451
## 18   45.35   64.55  1.14 0.3836604
## 19   40.81   63.68  1.47 0.2658518
## 20   29.65   36.03  1.11 0.2347251
## 21   37.07   46.18  2.41 0.1387512
## 22   43.54   53.99  3.86 0.1435462
## 23   48.30   59.64  4.47 0.1877425
## 24   67.55   93.43 13.02 0.2115034
## 25   27.62   31.76  0.58 0.4062370
## 26   34.85   40.87  1.94 0.1546763
## 27   42.35   50.92  3.78 0.1357390
## 28   48.53   59.08  4.51 0.1937344
## 29   54.70   68.64  5.12 0.2589905
## 30   57.50   73.24  3.98 0.3325961
## 31   66.40   89.07  5.02 0.3891421
plotSize(myFit2) # model fit to the data

It can be seen that the plots of the model fit for years 1994, 1999 and 2001 have superimposed ‘unrealistically high’ warnings. This is a subjective judgement warning made by the developer, and could be due to, for example, the F/M ratio exceeding 5, or the selectivity parameters being close to or above Linf.

Often these unrealistically high values may be caused by strong cohorts moving through the size compositions data. While this is not strongly apparent in our test data set, it is not uncommon - but is just about a perfect example of a data set that would violate LBSPR equilibrium assumptions. In such cases, users should be cautioned against using LBSPR or other similar equilibrium based approaches.

plotMat(myFit2) # specified maturity-at-length curve, and the estimated selectivity-at-length curve

plotEsts(myFit2) # Provides a plot of the estimated parameters

Here we can clearly see the unrealistically high selectivity values in 1994 and 2001, and unrealistically high (>5) F/M values in 1994, 1999 and 2001.

Interpreting results in a management context

The estimated values in gives you important management information.

The first is F/M, which is a proxy for overfishing where the assumption is that the fishing mortality at Maximum Sustainable Yield is assumed to be equal to natural mortality if recruitment is constant (Francis, 1974).

The other value is the SPR for the length data provided. This can be compared to the target SPR (in our case 40%) and the ratio SPR/SPRtarget gives whether overfishing is occurring or not.

Key references

LBSPR Methodology

Hordyk, A.R., Ono, K., Sainsbury, K.J., Loneragan, N. and J.D. Prince. 2015a. Some explorations of the life history ratios to describe length composition, spawning-per-recruit, and the spawning potential ratio. ICES Journal of Marine Science 72:204-216. https://doi.org/10.1093/icesjms/fst235.

Hordyk, A.R., Ono, K., Valencia, S.R., Loneragan, N.R. and J.D. Prince. 2015b. A novel length-based empirical estimation method of spawning potential ratio (SPR), and tests of its performance, for small-scale, data-poor fisheries. ICES Journal of Marine Science 72:217-231. https://doi.org/10.1093/icesjms/fsu004.

Hordyk, A.R., Loneragan, N.R. and J.D. Prince. 2015c. An evaluation of an iterative harvest strategy for data-poor fisheries using the length-based spawning potential ratio assessment methodology. Fisheries Research 171:20-32. https://doi.org/10.1016/j.fishres.2014.12.018.

Hordyk, A., Ono, K., Prince, J.D. and C.J. Walters. 2016. A simple length-structured model based on life history ratios and incorporating size-dependent selectivity: application to spawning potential ratios for data-poor stocks. Canadian Journal of Fisheries and Aquatic Sciences 73:1787-1799. https://doi.org/10.1139/cjfas-2015-0422.

Prince, J.D., Hordyk, A.R., Valencia, S.R., Loneragan, N.R. and K.J. Sainsbury. 2015. Revisiting the concept of Beverton-Holt life-history invariant with the aim of informing data-poor fisheries assessment. ICES Journal of Marine Science 72:194-203. https://doi.org/10.1093/icesjms/fsu011.

LBSPR Application

Prince, J.D., Victor, S., Kloulchad, V. and A.R. Hordyk. 2015. Length based SPR assessment of eleven Indo-Pacific coral reef fish populations in Palau. Fisheries Research 171:42-58. https://doi.org/10.1016/j.fishres.2015.06.008.

Other

Clark, W.G. 1993. The effect of recruitment variabilityon the choice of a target level of spawning biomass per recruit. Pages 233–246 in G. Kruse, D. M. Eggers, R. J. Marasco, C. Pautzke, and T. J. Quinn II, editors. Proceedings of the international symposium on management strategies for exploited fish populations. University of Alaska, Alaska Sea Grant Report 93-02., Fairbanks.

Clark, W.G. 2002. F35% revised ten years later. em>North American Journal Fisheries Management. https://doi.org/10.1577/1548-8675(2002)022<0251:FRTYL>2.0.CO;2.

Dorn, M.W. 2002. Advice on west coast rockfish harvest rates from Bayesian meta-analysis of stock–recruit relationships. North American Journal Fisheries Management 22:280–300. https://doi.org/10.1577/1548-8675(2002)022<0280:AOWCRH>2.0.CO;2

Francis, R.C. 1974. Relationship of fishing mortality to natural mortality at the level of maximum sustainable yield under the logistic stock production model. Journal of the Fisheries Research Board of Canada 31(9):1539–1542. https://doi.org/10.1139/f74-189.

Mace, P.M. and. M.P. Sissenwine. 1993. How much spawning per recruit is enough? Pages 101–118 in G. Kruse, D. M. Eggers, R. J. Marasco, C. Pautzke, and T. J. Quinn II, editors. Proceedings of the international symposium on management strategies for exploited fish populations. University of Alaska, Alaska Sea Grant Report 93-02., Fairbanks.

Myers, R.A., Barrowman, N.J., Hilborn, R. and D.G. Kehler. 2002. Inferring Bayesian priors with limited direct data: Applications to risk analysis. North American Journal Fisheries Management 22:351-264. https://doi.org/10.1577/1548-8675(2002)022<0351:IBPWLD>2.0.CO;2.