Catch curves (with selectivity)

Catch curve analysis is a method for estimating the total mortality of a stock (Z). The rate at which individuals die can be estimated from the slope ofa regression of the logarithms of relative numbers present in each age class on age. It can be used whenever there is one or more years of catch-at-age data (or at-length data, if they can be converted to age). The data can be fishery-dependent or -independent so long as they are representative of the population’s relative age/length structure. Fishing mortality (F) can be estimated as Z-M given an estimate of natural mortality (M) from another source (e.g., the literature, from a marine protected area, or from a tagging study).

The method is based on the equilibrium assumption (i.e., the total population receives constant recruitment and there is a constant mortality rate each year), which implies that all cohorts would be identical and the numbers in each age class would decline exponentially. Thus, if one assumes that fishing mortality and natural mortality have been constant, and one has the age-composition of the catch for a single year, if one plots the log-transformed numbers-at-age against age, then the slope of a linear regression through the data will provide an estimate of total mortality (Z).

Conventional catch-curve methods rely on the strong assumptions of constant fishing and natural mortality rates above some fully selected age that is usually estimated by visually inspecting a plot of log-catch at age vs age. Rather than just ignoring the earlier weakly selected ages below the (assumed) fully selected age, as in the classical catch-curve, it is possible to estimate selectivity parameters as well as the total mortality rate. However, this requires the use of an equilibrium based age-structured model rather than the log-transformed linear regression (see equations in Wayte and Klaer, 2010).

Catch curves with selectivity use a simple age-structured model to estimate selectivity from the age data. This provides an estimate of fully-selected fishing mortality rather than an average fishing mortality applied to all included age classes (as in the classic catch curve). This also circumvents the need to select a minimum age from which to fit the regression (to avoid the effects of selectivity), as is required for the classical catch curve.

As opposed to classic catch curves, catch curves with selectivity provide an estimate of fully-selected fishing mortality, rather than an average fishing mortality applied to all included age classes.

Although there are several packages that can be used to run catch curves with selectivity, the one by Dr Malcolm Haddon is supported here, because his version has been applied in several cases in Australia (refs needed).

See the Appendix for further description and equations, as supplied by Dr. Haddon and referencing Wayte and Klaer (2010) (Fisheries Research 106: 310-320).

Installation

This package is run from R as found on CRAN. R can be installed from 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/.

The supported version of Catch Curve Analysis is part of a set of stock assessment tools within Dr Haddon’s datalowSA package. See the link https://github.com/haddonm/datalowSA.

Therefore one needs to first install datalowSA from within R.

Set up checking/installing packages functions

# Ordinary checking/installing function
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)
}

# github checking/installing function
Github_CheckInstall.packages <- function(pkg,github_pkg,build_vignettes){
  new.pkg <- pkg[!(pkg %in% installed.packages()[, "Package"])]
  if (length(new.pkg)) 
    devtools::install_github(pkg,github_pkg,build_vignettes)
  sapply(pkg, require, character.only = TRUE)
}

Installing packages from GitHub directly.

Check and install “devtools”

Check_Install.packages("devtools")
## [1] "'devtools' has been installed already!"
## Loading required package: devtools
## Loading required package: usethis
## devtools 
##     TRUE

Check and install datalowSA from github

# Check if existed and if not, install the required package TMBhelper:
github_pkg = "haddonm/datalowSA"
Github_CheckInstall.packages("datalowSA",github_pkg=github_pkg,build_vignettes=TRUE)
## Loading required package: datalowSA
## datalowSA 
##      TRUE
#install_github("haddonm/datalowSA")

NB Once you have installed datalowSA, you don’t need to run the install.packages code chunk again. Luckily it will just give an error stating it won’t reinstall if your version is up to date.

Note that you may get an error, because the current version you have of Rcpp may preclude the installation of datalowSA. If this occurs, try reinstalling Rcpp from CRAN (Tools -> Install Packages).

library(datalowSA)

It should be noted that there are also good vignettes for this package by typing browseVignettes(package = “datalowSA”).

Reading in the data

To run the Catch Curve Analysis with Selectivity requires the following inputs:

  • Numbers-at-age * - for a single year, as a vector. This must be absolute numbers (as opposed to proportions-at-age) and a representative (across time) random sample.
  • A single natural mortality estimate *
  • The maximum age (within the available data) *
  • Age-at-50%-selectivity, a50 *
  • Delta (age-at-95% selectivity (a95) less the age-at-50% selectivity (a50), ie (a95-a50) *
  • An estimate of current fishing mortality (used as an initial value) *

(Note that delta feeds into the equation for selectivity-at-age:

Selectivity at age is modelled using (1+exp(-ln(19)(a - a50)/delta))-1 instead of (1+exp(-ln(19)(a -a50)/(a95-a50)))-1 to reduce the correlation between two estimated parameters (a50 and a95), and to avoid a50 values larger than a95 values. It is relatively easy to get a a95 value less than the a50 value, and this leads to a selectivity pattern that declines instead of increasing with age.

In general, it goes without saying that one requires a minimum number of ages in the data set to ensure that there are more data points than parameter values. The selectivity-included method in Dr. Haddon’s package estimates three parameters and it is a reasonable rule of thumb to have at least three data points per parameter estimated - so, at least nine ages here.

A simulated example dataset has been created for the whole toolbox. Let us now read in these data. We assume you have saved the data to the same directory as your working directory. 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.

The simulated catch-at-age data, TotalCatchAge, has years (1960-2001) in rows (2-43), with 25 age bins (from ages 1 to 25) in columns (B to Z). The first 11 years of this file have no data, so we will focus on year 1971 onwards.

Let us read in the data:

TotCatAge <- read.csv("TotalCatchAge.csv")
print(TotCatAge[12:42,])
##    YearsOut          A1           A2           A3          A4          A5
## 12     1971    1.910415    0.3925469    0.8258727   16.170486    1.877880
## 13     1972    5.840126    6.0808642    1.4250714    2.104864   33.731255
## 14     1973   62.394002   15.8097160   44.9978582    7.727642    9.657692
## 15     1974   24.236550  118.0074654   23.0015271   51.392508    8.639112
## 16     1975   55.103716   42.7568940  325.3436211   48.613787  102.768715
## 17     1976   12.188843  102.7895846   63.0607343  371.431144   55.189997
## 18     1977   55.643306   12.7038754  223.0947889  112.513335  527.059248
## 19     1978   76.643167   70.9197949   25.9093006  258.989382  109.240022
## 20     1979   81.047009  140.5462582  130.8785724   40.843807  371.659430
## 21     1980   26.807254  138.8902626  231.3314230  156.984856   37.997537
## 22     1981  113.402232   32.9076055  294.9740543  334.191571  210.562848
## 23     1982  183.714672  132.7821316   59.2440875  379.157799  387.005313
## 24     1983  148.349856  389.6158453  257.9485600   96.057167  504.854930
## 25     1984   85.834605  363.2626107  686.4943428  349.365257   95.529260
## 26     1985  133.084778   75.6188777  470.2082807  769.466339  298.351901
## 27     1986  101.715251  259.6285064  119.7158847  641.722041  785.347181
## 28     1987  286.446263  135.5081637  443.1192489  172.312856  668.293122
## 29     1988  308.013299  455.4089268  248.9493585  581.282370  167.764770
## 30     1989  596.622027  473.7542170  599.7932138  214.131427  384.245819
## 31     1990 1122.484109  373.2161377  623.8125995  511.509113  159.229529
## 32     1991  148.449770 1989.3673138  520.7658285  559.171268  354.863224
## 33     1992  201.162806  123.3779336 2675.6240388  462.689508  388.971831
## 34     1993   27.397110  329.6243051  144.2653213 2151.587635  235.998073
## 35     1994  443.581717   50.0309362  477.3122639  158.166694 1697.555787
## 36     1995 1265.225845   18.0096479   89.0893601  431.900777  100.770044
## 37     1996   44.238434 2373.0474340   30.0203166   78.552392  277.098045
## 38     1997   31.438088    8.8899254  845.7349405    7.027043   11.366910
## 39     1998   16.823783   35.4404509    9.8173442  658.112337    3.784587
## 40     1999   10.469020   23.9587229   41.2237479    8.143751  520.974613
## 41     2000   26.843251   16.8107275   46.2265698   64.155627    8.694273
## 42     2001   39.994189   49.1654268   20.8476087   54.812648   67.803887
##             A6          A7         A8         A9         A10          A11
## 12   0.4235131   1.0764674   5.620313   1.527880   3.0488564 2.074872e+00
## 13   3.6735822   0.6932627   1.827491   9.158322   2.2808116 4.245410e+00
## 14 156.4327838  14.2957686   2.746602   8.506518  39.8109582 1.177706e+01
## 15   8.5749914 133.8152545  11.737072   3.696920   8.0920327 2.796508e+01
## 16  16.0962080  11.2722323 215.224237  23.309309   3.0988263 9.996217e+00
## 17  86.1130907  13.8538426  13.588940 160.812182  15.5758392 4.213590e+00
## 18  67.2522226 103.2832831  20.280873  10.532464 188.4009945 1.636358e+01
## 19 475.3970046  51.7767864  73.055668  10.190807  14.6019250 1.449906e+02
## 20 131.0130449 540.0044436  58.489899  80.579764  14.4076450 1.141600e+01
## 21 324.3874383 121.9559262 455.457120  46.251419  81.7177133 1.382336e+01
## 22  39.6009276 378.1201915 130.873968 453.578522  49.0770492 6.207193e+01
## 23 207.4801095  39.2712029 300.558800  87.149002 357.1745328 3.692002e+01
## 24 435.9721487 206.2606835  40.758918 270.312910  81.1141696 3.124874e+02
## 25 454.1988904 346.5794458 164.564215  31.421909 190.2727395 5.145767e+01
## 26  70.7423726 289.9335413 214.307032  91.564524  15.5194472 1.275537e+02
## 27 239.0346777  52.9139342 217.612957 143.570593  46.8706752 1.121730e+01
## 28 657.9779811 174.7945636  39.602581 132.139531  86.9198092 3.954351e+01
## 29 500.6156585 405.2064002 117.069545  16.773959  76.8019004 3.785172e+01
## 30  89.7955500 201.4489260 170.623425  38.231450   8.0747136 2.196131e+01
## 31 192.4781307  34.0222976  71.736309  40.816773   9.8335650 1.415969e+00
## 32  81.9139204  83.1418089  11.063863  23.858871  14.1746997 9.953650e-01
## 33 156.6798065  23.9010080  30.579305   3.196548   4.0268565 2.220936e+00
## 34 198.3413233  60.9602678  13.961915   5.095099   0.5160038 1.962090e+00
## 35 167.8767647  81.8787208  32.601248   5.236099   2.4601533 7.193254e-02
## 36 888.0046035  57.3773141  29.131254  15.749575   0.3869628 3.428750e-01
## 37  46.2954214 342.7527245  21.530982  13.586527   3.1047126 3.315267e-02
## 38  39.1589728   5.9501188  32.787798   1.514699   0.9710801 2.169095e-03
## 39   7.2952123  21.4419043   3.213071  14.830530   1.5364099 5.703590e-01
## 40   2.8247573   5.2037954  13.468349   1.851832   8.3255664 6.372236e-01
## 41 630.0968344   1.9245927   4.117840  13.754784   1.8336675 1.018619e+01
## 42  11.0778164 500.8286221   2.348458   5.200383   9.1653060 1.822362e+00
##             A12          A13          A14          A15          A16         A17
## 12 1.444049e+00 1.131696e+00 1.036129e+00 7.738659e-01  0.573972127  0.56727849
## 13 3.819904e+00 2.099353e+00 1.961948e+00 1.375408e+00  1.241416222  0.84977139
## 14 1.635968e+01 1.433114e+01 9.040685e+00 7.690224e+00  5.649773644  4.62796753
## 15 7.694873e+00 1.353481e+01 8.599677e+00 6.845094e+00  5.735573364  4.20218104
## 16 4.684429e+01 1.393139e+01 2.032851e+01 1.870078e+01 10.442120469  6.52826629
## 17 9.760003e+00 3.420136e+01 9.591182e+00 1.105961e+01 12.360166365  6.80116649
## 18 2.459466e+00 6.483654e+00 3.863917e+01 1.001368e+01 18.137117691 18.87582202
## 19 1.404280e+01 2.716079e+00 8.449055e+00 2.791701e+01  6.869512568 12.90826157
## 20 1.398043e+02 1.176627e+01 2.629969e+00 5.158894e+00 34.898452637 10.92050639
## 21 7.572340e+00 1.075001e+02 9.273752e+00 2.249562e+00  5.755758139 19.71812493
## 22 6.123684e+00 5.937859e+00 1.016002e+02 8.963431e+00  1.425079025  5.05897258
## 23 4.262454e+01 7.114140e+00 3.271413e+00 8.007298e+01  8.059648028  1.45956723
## 24 3.341098e+01 5.325094e+01 6.518383e+00 3.457130e+00 67.419743623  4.62477993
## 25 2.203049e+02 2.303407e+01 2.325284e+01 4.014672e+00  4.365335251 53.00089502
## 26 3.228934e+01 1.260496e+02 1.141339e+01 1.347403e+01  3.276615990  0.83420585
## 27 6.315319e+01 2.610839e+01 8.053372e+01 7.914731e+00 10.810479916  0.88292178
## 28 2.966183e+00 4.114913e+01 1.193642e+01 4.375533e+01  2.168250729  4.60665875
## 29 1.225767e+01 2.186655e+00 1.274821e+01 5.100118e+00 20.649468952  1.64731808
## 30 1.353283e+01 3.374532e+00 4.749975e-01 3.941835e+00  1.428536084  3.24862709
## 31 6.161638e+00 2.600022e+00 2.052710e+00 6.968363e-02  0.647079702  0.03371932
## 32 1.603187e-01 9.486003e-01 5.317343e-02 3.653987e-01  0.233024542  0.85584157
## 33 1.163430e+00 0.000000e+00 2.860974e-01 2.753646e-02  0.000000000  0.00000000
## 34 5.221866e-02 3.854199e-02 0.000000e+00 0.000000e+00  0.000000000  0.00000000
## 35 9.900001e-03 7.005426e-01 2.928362e-04 8.694902e-02  0.000000000  0.00000000
## 36 6.112787e-01 3.601882e-03 1.368553e-03 5.492618e-03  0.000000000  0.00000000
## 37 5.974887e-02 0.000000e+00 0.000000e+00 1.087493e-04  0.000000000  0.00000000
## 38 2.553828e-04 0.000000e+00 0.000000e+00 0.000000e+00  0.000000000  0.00000000
## 39 4.262780e-01 0.000000e+00 3.206832e-05 3.062924e-02  0.000000000  0.00000000
## 40 2.337456e-01 2.467208e-02 1.322278e-04 7.243238e-05  0.000000000  0.00000000
## 41 1.295448e+00 9.046427e-02 9.394603e-03 0.000000e+00  0.004335718  0.00000000
## 42 5.584710e+00 6.179673e-01 2.976959e-02 1.099712e-05  0.005835920  0.00000000
##            A18          A19          A20          A21          A22          A23
## 12  0.46707001 3.595825e-01 3.090172e-01 2.259975e-01 0.1676614955 1.446701e-01
## 13  0.81118788 7.796925e-01 5.400465e-01 3.191857e-01 0.3377839656 3.420404e-01
## 14  3.54970605 2.728570e+00 2.763659e+00 2.056554e+00 1.6854164047 1.310860e+00
## 15  3.98849818 3.057752e+00 2.824627e+00 2.501577e+00 1.1099254146 1.714971e+00
## 16  6.19360784 5.240554e+00 6.359342e+00 4.970663e+00 2.7921090204 2.941452e+00
## 17  5.56515872 4.332396e+00 4.527929e+00 3.956677e+00 2.5235344863 2.463087e+00
## 18  8.54163295 5.191591e+00 6.337373e+00 6.387457e+00 2.7464261142 2.100510e+00
## 19  8.67620315 6.959187e+00 5.840007e+00 4.688535e+00 3.6037951718 3.430368e+00
## 20  8.44497160 1.218222e+01 6.318749e+00 5.731237e+00 3.3787361919 3.267615e+00
## 21  7.17024165 6.457755e+00 7.691318e+00 7.184870e+00 3.1594607077 3.442483e+00
## 22 17.90600506 5.250411e+00 9.833281e+00 1.443184e+01 3.7422225763 2.533574e+00
## 23  1.50802738 1.253483e+01 4.699851e+00 6.034106e+00 5.4671850727 1.664063e+00
## 24  1.28162311 4.080661e+00 1.275567e+01 6.221721e+00 3.3483318441 4.909587e+00
## 25  3.30103442 2.984081e-01 3.401570e+00 8.858118e+00 3.1794387167 5.218498e+00
## 26 22.72123598 1.731936e+00 8.869921e-02 1.342427e+00 6.0933784180 1.239761e+00
## 27  0.70708219 1.851978e+01 8.016094e-01 4.756846e-02 0.6436520564 1.835686e+00
## 28  0.46554735 1.294484e-02 8.188379e+00 9.172146e-02 0.0004973671 6.203893e-01
## 29  1.17337409 3.994768e-02 1.588208e-01 1.670424e+00 0.1651683882 0.000000e+00
## 30  0.35092115 4.762157e-01 0.000000e+00 1.782347e-04 0.4652783693 1.126879e-04
## 31  0.39929430 4.009703e-03 2.264478e-01 0.000000e+00 0.0000000000 7.112896e-03
## 32  0.00000000 1.384126e-04 5.921586e-04 6.979977e-03 0.0000000000 0.000000e+00
## 33  0.05235793 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 34  0.00000000 9.666143e-05 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 35  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 36  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 37  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 38  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 39  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 7.157566e-05
## 40  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 41  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
## 42  0.00000000 0.000000e+00 0.000000e+00 0.000000e+00 0.0000000000 0.000000e+00
##           A24         A25
## 12 0.08024316  0.45846617
## 13 0.35627317  1.11742859
## 14 0.80927973  4.71588408
## 15 0.72236562  3.50056611
## 16 2.71981575  9.06333014
## 17 2.04740905  5.55253487
## 18 2.81083035 10.01730059
## 19 1.98934474 10.12539818
## 20 1.10508659  9.83712957
## 21 2.45434518 10.49553551
## 22 2.32453925 13.20799849
## 23 2.04154336 10.50043192
## 24 3.04258550 17.19525948
## 25 3.51053076  9.28874743
## 26 3.92435221  6.35022780
## 27 0.60345904  9.05873048
## 28 1.25969016  3.92123258
## 29 0.16992485  0.92499674
## 30 0.00000000  0.86386005
## 31 0.00000000  0.04374913
## 32 0.00000000  0.00000000
## 33 0.00000000  0.00000000
## 34 0.00000000  0.00000000
## 35 0.00000000  0.00000000
## 36 0.00000000  0.00000000
## 37 0.00000000  0.00000000
## 38 0.00000000  0.00000000
## 39 0.00000000  0.00000000
## 40 0.00000000  0.00000000
## 41 0.00000000  0.00000000
## 42 0.00000000  0.00000000

On inspection, there is no catch within many older ages from 1988. On the other hand, the earlier years also appear to have very few counts relative to those in later years.

We can visualise the age composition for specified years using the function compyear. In our case, we will focus on 1980-1995 (rows 21 to 36) inclusive. (Note that, when plotting, we will have to remove the year column from the data).

naa <- TotCatAge
yrs <- 21:36
naa <- naa[yrs,2:26] # We must remove the first column containing the years
ages <- 1:25 
#plotprep() 
compyear(naa,1980:1995,ages,plots=c(4,4),freq=TRUE,border=3) 

As per the advice in Dr. Haddon’s vignettes: when our simulated example age data are considered using these frequency plots, it is clear that the notion of a exponential decline is only approximated by these data. Such deviations from the theoretical expectation derive from the fact that real populations are not in equilibrium and recruitment is not constant. So the dynamics appear to be realistic although the assumptions about equilibrium do not appear to be valid. Despite this we can use these data to illustrate the use of catch curves.

The sample sizes from each year are printed on the above plot, but can also be found by using the following.

sample_sizes <- apply(TotCatAge[,2:26],1,sum)
sample_sizes
##  [1]    0.0000    0.0000    0.0000    0.0000    0.0000    0.0000    0.0000
##  [8]    0.0000    0.0000    0.0000    0.0000   42.6888   87.0125  451.4760
## [15]  485.1910 1010.6400 1013.5600 1475.8700 1429.9300 1856.3300 1835.7300
## [22] 2297.7000 2357.5100 2965.2500 3184.0100 2787.1800 2840.9700 2957.8000
## [29] 2974.6300 2826.8400 3152.8000 3790.3900 4073.9600 3169.8000 3117.5700
## [36] 2896.6100 3230.3200  984.8420  773.3230  637.3400  826.0440  769.3050

We can see a strong year cohort progressing through from 1990 onwards, which violates the assumption of constant recruitment. Despite the large sample sizes in the 1990s, let us instead focus our analysis on years 1986 and 1988 in the first instance. These are years with high sample sizes (2840 and 2974, respectively), but without visual evidence of strong recruitment pulses.

Our dataset was created using Stock Synthesis based on an age-at-maturity of roughly 3-4 years and a rate of natural mortality of 0.2 yr -1 for females, and 0.25 yr -1 for males.

Regarding the input parameters:

  • Natural mortality estimate = 0.225 yr -1 (average of males and females)
  • The maximum age in the sample data = 25.
  • Age-at-50%-selectivity = 4.45 years (average of males and females: 3.97 for females, 4.93 for males)
  • Age-at-95%-selectivity = 19.07 years (average of males and females: 13.14 for females, 19.07 for males)
  • Delta = 19.07-4.45 = 14.62 years
  • An estimate of current fishing mortality (used as an initial value) =

    ???STILL NEEDED

Fitting the model to our data

selectCC is the function that sets up and runs the catch curve (with selectivity) model and simulations.

The main input data expected is a vector of numbers-at-age for a single year. As stated above, we will try years 1986 and 1988 (rows 27 and 29 of TotCatAge), separately.

We will assign the same names to the numbers-at-age vector and the other inputs as per datalowSA.

In creating our catch-at-age vector, there are a few potential trip-points. The first is that the data (numbers-at-age, numaa) MUST be a vector, not a list. The second is that the rownames (or names) of the vector should equate to a character list of the actual ages.

numaa <- TotCatAge[27,2:26]
names(numaa) <- as.character(1:25)
numaa <- unlist(numaa)
print(numaa)
##            1            2            3            4            5            6 
## 101.71525140 259.62850639 119.71588468 641.72204084 785.34718070 239.03467769 
##            7            8            9           10           11           12 
##  52.91393420 217.61295704 143.57059311  46.87067516  11.21729739  63.15319038 
##           13           14           15           16           17           18 
##  26.10838594  80.53372129   7.91473118  10.81047992   0.88292178   0.70708219 
##           19           20           21           22           23           24 
##  18.51978373   0.80160937   0.04756846   0.64365206   1.83568558   0.60345904 
##           25 
##   9.05873048

The model needs the natural mortality and maximum sampled age as a globals list which we can call glb, copying the example in datalowSA. Let’s create the list for glb (which we could of course have instead read in).

NatM <- 0.225
maxage <- 25
glb <- list(M = NatM,maxage = maxage)
print(glb)
## $M
## [1] 0.225
## 
## $maxage
## [1] 25

Finally, the parameters relating to selectivity:

pars <- c(A50=4.45,delta=14.62,fcur=0.6)
print(pars)
##   A50 delta  fcur 
##  4.45 14.62  0.60

Running the model

Now that we have the data and inputs set up, the next is to run the model itself, using selectCC.

library(datalowSA)
out <- selectCC(glb$M,glb$maxage,numaa,pars,plot=TRUE)

This shows the fitted catch curve.

We can examine the model output to check the fitted selectivity parameters (age-at-50%-maturity, A50, delta, and current fishing mortality, fcur) and predicted catch-at-age:

print(out)

We can see that, for 1986, this results in a fully selected fishing mortality estimate of 0.1342896 yr-1 (and a logistically-increasing selectivity-at-age, as indicated in the resultant plot).

Using instead the catch-at-age data from 1988:

numaa <- TotCatAge[29,2:26]
names(numaa) <- as.character(1:25)
numaa <- unlist(numaa)
print(numaa)
##            1            2            3            4            5            6 
## 308.01329936 455.40892679 248.94935853 581.28236991 167.76477039 500.61565845 
##            7            8            9           10           11           12 
## 405.20640019 117.06954536  16.77395913  76.80190036  37.85171579  12.25767033 
##           13           14           15           16           17           18 
##   2.18665527  12.74820878   5.10011753  20.64946895   1.64731808   1.17337409 
##           19           20           21           22           23           24 
##   0.03994768   0.15882083   1.67042422   0.16516839   0.00000000   0.16992485 
##           25 
##   0.92499674
out <- selectCC(glb$M,glb$maxage,numaa,pars,plot=TRUE)

print(out)
## $best
## $best$par
##       A50     delta      fcur 
## 4.0913906 3.6463649 0.2494271 
## 
## $best$value
## [1] 7061.287
## 
## $best$counts
## function gradient 
##      244       NA 
## 
## $best$convergence
## [1] 0
## 
## $best$message
## NULL
## 
## 
## $result
##          counts   predcounts    logcount logpredcount         pCAA     predpCAA
## 1  308.01329936 196.54044710  5.73014296    5.2808682 1.035468e-01 6.607223e-02
## 2  455.40892679 306.30674909  6.12119576    5.7245870 1.530977e-01 1.029731e-01
## 3  248.94935853 420.17552380  5.51724950    6.0406725 8.369087e-02 1.412530e-01
## 4  581.28236991 478.40472749  6.36523665    6.1704571 1.954133e-01 1.608283e-01
## 5  167.76477039 442.66904428  5.12256282    6.0928224 5.639853e-02 1.488148e-01
## 6  500.61565845 344.91370445  6.21583866    5.8432943 1.682951e-01 1.159518e-01
## 7  405.20640019 240.59619239  6.00439657    5.4831200 1.362208e-01 8.088273e-02
## 8  117.06954536 158.07088457  4.76276816    5.0630436 3.935600e-02 5.313968e-02
## 9   16.77395913 100.84666966  2.81982763    4.6136012 5.639007e-03 3.390226e-02
## 10  76.80190036  63.46394799  4.34122938    4.1504720 2.581898e-02 2.133507e-02
## 11  37.85171579  39.69062322  3.63367631    3.6811150 1.272485e-02 1.334305e-02
## 12  12.25767033  24.75321317  2.50615189    3.2089553 4.120738e-03 8.321443e-03
## 13   2.18665527  15.41807736  0.78237310    2.7355407 7.351016e-04 5.183192e-03
## 14  12.74820878   9.59810183  2.54539077    2.2615654 4.285645e-03 3.226654e-03
## 15   5.10011753   5.97354035  1.62926358    1.7873398 1.714538e-03 2.008162e-03
## 16  20.64946895   3.71731841  3.02768960    1.3130026 6.941861e-03 1.249674e-03
## 17   1.64731808   2.31316225  0.49914856    0.8386155 5.537892e-04 7.776302e-04
## 18   1.17337409   1.43937111  0.15988344    0.3642063 3.944605e-04 4.838824e-04
## 19   0.03994768   0.89564348 -3.22018458   -0.1102128 1.342946e-05 3.010941e-04
## 20   0.15882083   0.55730846 -1.83997855   -0.5846364 5.339179e-05 1.873539e-04
## 21   1.67042422   0.34678096  0.51307762   -1.0590619 5.615570e-04 1.165795e-04
## 22   0.16516839   0.21578163 -1.80078979   -1.5334883 5.552569e-05 7.254066e-05
## 23   0.00000000   0.13426831          NA           NA 0.000000e+00 4.513782e-05
## 24   0.16992485   0.08354732 -1.77239897   -2.4823421 5.712470e-05 2.808663e-05
## 25   0.92499674   0.13761861 -0.07796507   -1.9832691 3.109619e-04 4.626411e-05
##    Selectivity        FaA
## 1   0.07611816 0.01898593
## 2   0.15593405 0.03889418
## 3   0.29290894 0.07305943
## 4   0.48155891 0.12011385
## 5   0.67561712 0.16851723
## 6   0.82363936 0.20543799
## 7   0.91283097 0.22768479
## 8   0.95915238 0.23923861
## 9   0.98136131 0.24477812
## 10  0.99160094 0.24733216
## 11  0.99623675 0.24848846
## 12  0.99831819 0.24900763
## 13  0.99924926 0.24923986
## 14  0.99966505 0.24934357
## 15  0.99985059 0.24938985
## 16  0.99993336 0.24941049
## 17  0.99997028 0.24941970
## 18  0.99998675 0.24942381
## 19  0.99999409 0.24942564
## 20  0.99999736 0.24942646
## 21  0.99999882 0.24942682
## 22  0.99999948 0.24942698
## 23  0.99999977 0.24942706
## 24  0.99999990 0.24942709
## 25  0.99999995 0.24942710

For 1988, fully selected fishing mortality is estimated to be 0.2494271 yr-1 (again with a logistically-increasing selectivity-at-age, indicated in the resultant plot).

(Note: using 1990’s data (row 31) results in a fishing mortality estimate of 0.1949766 yr-1 for all ages. The fitted proportion-at-age plot shows a negative exponential decline, which equates to all age classes being fully selected. This is the result of the strong cohort in the earliest age class. Using 1992’s data (row 33), we see an increase, a peak and then a decline in proportion-at-age (reflecting the strong cohort having since aged by 2 years), which yielded a more logistic increase in selectivity with age. The fully-selected fishing mortality was estimated to be very high, at 0.952 yr-1, again as a result of this strong year cohort).

Other years can be tested, as desired. It is important to reiterate that the data do not result in a constant fishing mortality over time.

Interpreting results in a management context

Catch curve analysis results in an estimate of total mortality, and, assuming a fixed natural mortality, fishing mortality. It does not provide an estimate of stock status per se, except where consistently low fishing mortality is used as a proxy for sustainability (although this would require many years of age-composition data). The estimates of fishing mortality could be compared to a reference point such as the fishing mortality at maximum sustainable yield, if this were available.

As only a single year of numbers-at-age is required, catch curve analysis suits fisheries which have occasional age-structured data, and whose catch-at-age distribution is not expected to vary strongly through time. As such, catch-curve analysis is frequently applied to data-poor fisheries, where a variety of disparate information sources (each of low quality for stock assessment) may be available. Future work could formalize a process by which catch curves are combined with expert opinion (i.e., interviews with fishers) and other data-poor analyses (e.g., changes in mean length) into a quantitative estimate of stock status (Thorson and Prager 2011).

Key references

Catch Curve Methodology

Chapman, D.G. and D.S. Robson. 1960. The analysis of a catch curve. Biometrics 16:354-368.

Dunn, A., Francis, R.I.C.C. and I.J. Doonan. 2002. Comparison of the Chapman-Robson and regression estimators of Z from catch-curve data when non-sampling stochastic error is present. Fisheries Research 59:149-159. https://doi.org/10.1016/S0165-7836(01)00407-6

Gulland, J.A. 1971. The fish resources of the ocean. West Byfleet, UK: Fishing News Books.

Smith, M.W., Then, A.Y., Wor, C., Ralph, G., Pollock, K.H. and J.M. Hoenig. 2012. Recommendations for catch-curve analysis. North American Journal of Fisheries Management 32:956-967. http://dx.doi.org/10.1080/02755947.2012.711270

Thorson, J.T. and M.H. Prager. 2011. Better catch curves: Incorporating age-specific natural mortality and logistic selectivity. Transactions of the American Fisheries Society 140:356-366. http://dx.doi.org/10.1080/00028487.2011.557016

Wayte, S.E. and N.L. Klaer. 2010. An effective harvest strategy using improved catch-curves. Fisheries Research 106:310-320. https://doi.org/10.1016/j.fishres.2010.08.012

Catch Curve Evaluation

Allen, M.S. 1997. Effects of variable recruitment on catch-curve analysis for crappie populations. North American Journal of Fisheries Management 17(1):202-205. https://doi.org/10.1577/1548-8675(1997)017<0202:EOVROC>2.3.CO;2

Griffiths, S.P. 2010. Stock assessment and efficacy of size limits on longtail tuna (Thunnus tonggol) caught in Australian waters. Fisheries Research 102(3):248-257. https://doi.org/10.1016/j.fishres.2009.12.004

Oyarzún, C., Cortés, N. and E. Leal. 2013. Age, growth and mortality of southern rays bream Brama australis (Bramidae) off the southeastern Pacific coast. Revista de biología marina y oceanografía 48(3). http://dx.doi.org/10.4067/S0718-19572013000300014

See references in http://derekogle.com/fishR/examples/oldFishRVignettes/CatchCurve.pdf.

Appendix: Catch Curve Equations

Catch-Curves with Selectivity

(reproduced with permission from Dr. Haddon’s vignettes)

A simple model of the expected explonential decline of numbers-at-age is needed:

\[\begin{matrix} {{N}_{0}}= & 1 & a=0 \\ {{N}_{a}}= & {{N}_{a-1}}\exp \left( -{{s}_{a-1}}F-M \right) & 0<a<{{a}_{\max }} \\ {{N}_{a}}= & \frac{{{N}_{a-1}}\exp \left( -{{s}_{a-1}}F-M \right)}{\left( 1-\exp \left( -{{s}_{a}}F-M \right) \right)} & a={{a}_{\max }} \\ \end{matrix} \]

where \(N_a\) is the numbers at age \(a\), \(s_a\) is the selectivity at age \(a\), \(F\) is the fully- selected average fishing mortality over the years for which agedtaa are available, \(M\) is the assumed rate of natural mortality, and \(a_{max}\) is the maximum age. Selectivity is assumed to be asymptotic (so this method may not be suitable for data collected using gillnets or hook gear):

\[s_a = \frac{1}{1+exp(-\ln(19)(a-a_{50})/\delta)}\] where \(a\) is age, \(a_{50}\) is the age-at-50%-selectivity and \(\delta\) is equal to \(a_{95} - a_{50}\).

The use of 19 is something of an enigma: if you do not use 19, the model does not predict the 95% of the selectivity, but it is not a simple relationship.

It should also be noted that the selectivity curve used does not need to be the particular version of the logistic that is used herein. Instead one can use other logistic equations, and as they do not use the rather arcane -Ln(19) they are sometimes easier to understand

Three parameters are estimated, both the selectivity paramters, \(a_{50}\) and \(\delta\), and the average fully-selected \(F\) value. \(\delta\) is used instead of \(a_{95}\) as this reduces the correlation between the two selectivity paarameters and helps avoid \(a_{50}\) values larger than the \(a_{95}\) values. These parameters can be estimated through minimizing the multinomial negative log-likelihood:

\[-veLL=-\sum\limits_{y}\sum\limits_{a}{{O}_{y,a}}log\left(\frac{{\hat{N}_{a}}}{\sum{{\hat{N}_{a}}}} \right)\]

where \(O_{y,a}\) is the observed proportion of numbers-at-age \(a\) in the catch during year \(y\), and \({\hat{N}_{a}}\) is the expected or predicted proportion of numbers-at-age \(a\) in the catch during year \(y\). The predicted numbers-at-age are calculated using the Baranov catch equation:

\[{\hat{N}_a}=\frac{{N_a}{s_a}F}{(s_aF+M)}\left(1-e^{-(s_aF+M)} \right )\]