> library(FSA) # catch.curve
> library(FSA) # catch.curve
The data for this box is very simple — number of fish in ten age-classes — and could have been entered directly into R as shown in the Box6.2 vignette. However, I have entered these data into an external tab-delimited text file (same as that used in Box6.3 vignette) which is read into R as follows:
> setwd("c://aaaWork//web//fishR//bookex//AIFFD//Box6_4")
> d <- read.table("box6_3.txt",header=TRUE)
> d
age catch
1 1 90
2 2 164
3 3 162
4 4 110
5 5 55
6 6 41
7 7 20
8 8 14
9 9 7
10 10 5
The catch-curve methods illustrated below require the natural logarithm of catches. This variable is added to the data frame with,
> d$logcatch <- log(d$catch)
Linear regressions in R are performed with lm() where the first argument is a formula of the form y~x and the second (data=) argument is the data frame in which the variables can be found. An optional third (subset=) argument can be used to efficiently create a subset of data frame. This argument is particularly important in catch-curve analyses because it can be used to limit the analysis to only those ages that appear to be fully recruited to the gear. In this example, the authors have determined that all fish age-3 and older are fully recruited and should be used to estimate the mortality rate. The "catch-curve" regression is then fit with,
> cc1 <- lm(logcatch ~ age, data = d, subset = age >= 3)
The ANOVA table and parameter estimates are extracted with,
> anova(cc1)
Analysis of Variance Table
Response: logcatch
Df Sum Sq Mean Sq F value Pr(>F)
age 1 10.9766 10.9766 1072.5 5.391e-08
Residuals 6 0.0614 0.0102
> summary(cc1)
Call:
lm(formula = logcatch ~ age, data = d, subset = age >= 3)
Residuals:
Min 1Q Median 3Q Max
-0.11343 -0.08876 0.01113 0.07263 0.12057
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.66033 0.10758 61.91 1.19e-09
age -0.51122 0.01561 -32.75 5.39e-08
Residual standard error: 0.1012 on 6 degrees of freedom
Multiple R-squared: 0.9944, Adjusted R-squared: 0.9935
F-statistic: 1073 on 1 and 6 DF, p-value: 5.391e-08
with 95% default CIs for the parameters extracted with,
> confint(cc1)
2.5 % 97.5 %
(Intercept) 6.3970827 6.9235800
age -0.5494179 -0.4730256
As noted in Box 6.4, the estimated of the instantaneous mortality rate (Z) is obtained from the slope estimated. Thus, Z (and 95% CI) can be isolated with, [keep.source=TRUE] Z ← -coef(cc1)[2] # 2nd coefficient of model results Z ZCI ← -confint(cc1)[2,] # 2nd row of confint results ZCI
These results can be used to more effectively find estimates for the annual mortality rate (A),
> A <- 100 * (1 - exp(-Z))
> A
age
40.03373
> ACI <- 100 * (1 - exp(-ZCI))
> ACI
2.5 % 97.5 %
42.39725 37.57322
As noted in Box 6.4 the predictions from the unweighted regression fit above can be stored and used as weights in the catch-curve analysis to minimize the effect of few fish in the older age classes on the regression results. Making these predictions is simple but using them is slightly complicated because of the subsetting required for the catch-curve regression. The easiest way to handle this slight complexity is to predict log catches for all ages, attach these to the additional data frame, and then use the same subsetting routine. Predictions are made with predict() which requires the fitted model as the first argument and a data frame that has a variable that is named exactly as the variable in the model as the second argument. Thus, the predictions for all ages in the origianl data frame is made and appended to the data frame with.
> d$W <- predict(cc1, d) > d age catch logcatch W 1 1 90 4.499810 6.149110 2 2 164 5.099866 5.637888 3 3 162 5.087596 5.126666 4 4 110 4.700480 4.615444 5 5 55 4.007333 4.104223 6 6 41 3.713572 3.593001 7 7 20 2.995732 3.081779 8 8 14 2.639057 2.570557 9 9 7 1.945910 2.059336 10 10 5 1.609438 1.548114
The second weighted regression is fit and the results extracted with,
> cc2 <- lm(logcatch ~ age, data = d, subset = age >= 3, weights = W)
> anova(cc1)
Analysis of Variance Table
Response: logcatch
Df Sum Sq Mean Sq F value Pr(>F)
age 1 10.9766 10.9766 1072.5 5.391e-08
Residuals 6 0.0614 0.0102
> summary(cc2)
Call:
lm(formula = logcatch ~ age, data = d, subset = age >= 3, weights = W)
Residuals:
Min 1Q Median 3Q Max
-0.196521 -0.153502 -0.006143 0.128350 0.228641
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.66128 0.10002 66.60 7.71e-10
age -0.51139 0.01643 -31.12 7.31e-08
Residual standard error: 0.1822 on 6 degrees of freedom
Multiple R-squared: 0.9938, Adjusted R-squared: 0.9928
F-statistic: 968.3 on 1 and 6 DF, p-value: 7.314e-08
> Z <- -coef(cc2)[2]
> Z
age
0.5113879
> ZCI <- -confint(cc2)[2, ]
> ZCI
2.5 % 97.5 %
0.5515999 0.4711759
The methods described above are implemented with catch.curve() from the FSA package. This function requires the vector of ages as the first argument, the vector of catches as the second argument, and a third argument that indicates which ages represent fully recruited fish. The ages of fully recruited fish can be identified with a full range (i.e., 3:10) or a concatenated list of ages (i.e., c(3,4,5,6,7,8,9,10)). In addition, the weighted regression can be performed by including use.weights=TRUE in catch.curve(). The results of this constructor function should be saved to an object so that specific information can be extracted from it. For example,
> cc3 <- catch.curve(d$age, d$catch, 3:10)
> anova(cc3)
Analysis of Variance Table
Response: log.catch.e
Df Sum Sq Mean Sq F value Pr(>F)
age.e 1 10.9766 10.9766 1072.5 5.391e-08
Residuals 6 0.0614 0.0102
> summary(cc3)
Estimate Std. Error t value Pr(>|t|)
Z 0.5112218 0.01560993 32.74978 5.3913e-08
A 40.0237632 NA NA NA
> confint(cc3)
95% LCI 95% UCI
Z 0.4730256 0.5494179
A 37.6885900 42.2714237
In addition, a summary graphic is produced with,
> plot(cc3)
The weighted regression results are,
> cc4 <- catch.curve(d$age, d$catch, 3:10, use.weights = TRUE)
> summary(cc4)
Estimate Std. Error t value Pr(>|t|)
Z 0.5113879 0.01643377 31.11812 7.314447e-08
A 40.0337256 NA NA NA
> confint(cc4)
95% LCI 95% UCI
Z 0.4711759 0.5515999
A 37.5732237 42.3972487