The working directory is set, the data file (i.e., box5_6.txt) is read, and the structure is observed with,
> setwd("c://aaaWork//web//fishR//bookex//AIFFD//Box5_6")
> d <- read.table("box5_6.txt", header = TRUE)
> str(d)
'data.frame': 17 obs. of 5 variables:
$ markd : Factor w/ 16 levels "3/30/98","5/10/92",..: 6 12 15 15 13 9 11 1 5 3 ...
$ clmark : num 70.3 60.5 65.6 61.2 76.9 64.4 97.4 60.9 62.5 67 ...
$ recapd : Factor w/ 16 levels "11/15/95","12/24/94",..: 9 13 14 10 9 11 4 7 5 12 ...
$ clrecap : num 76 64.7 67.9 65.2 79.1 67.5 97.9 67.7 69.1 69.5 ...
$ timeoutd: int 1091 1090 714 1056 715 697 335 1544 1087 427 ...
As noted in Box 5.6, the timeoutd variable is the time-at-large in days and should be converted to time-at-large in years with,
> d$timeouty <- d$timeoutd/365
As noted in the Box 5.4 vignette the von Bertalanffy model, even with the Fabens method, is a non-linear model that will require initial values for the model fitting. Unlike with the traditional von Bertalanffy model illustrated in the Box 5.4 vignette graphic methods for finding the starting values have not been developed for Fabens' method. Fortunately, finding starting values is fairly straight-forward. A simple starting value for the asymptotic mean length is the maximum length in the data frame,
> svLinf <- max(d$clrecap)
> svLinf
[1] 97.9
A reasonable starting value for the Brody growth coefficient is, as described here, the average (from all fish in the sample) instantaneous growth rate between the time at marking and the time a recapture. This average is computed with,
> svK <- with(d, mean((log(clrecap) - log(clmark))/timeouty))
> svK
[1] 0.05174481
These intial values are then entered into a list as shown for the traditional von Bertalanffy model,
> Fsv <- list(Linf = svLinf, K = svK)
The Fabens method von Bertalanffy model can be placed into a formula object as follows,
> Fvbmdl <- clrecap ~ clmark + (Linf - clmark) * (1 - exp(-K * timeouty))
where you must make sure that clrecap and clmark are the exact names of the variables in your data frame and Linf and K are the exact names of the parameters that you stored in your starting values list.
The non-linear model fitting procedure in R is implemented with nls(), which requires the model formula, the list of starting values, and the data frame as arguments. The results of the model fit should be stored into an object. The Febens method von Bertalanffy model is fit to these data with,
> tvb <- nls(Fvbmdl, start = Fsv, data = d)
The coefficient estimates and the correlations among coefficient estimates are extracted with,
> summary(tvb, correlation = TRUE)
Formula: clrecap ~ clmark + (Linf - clmark) * (1 - exp(-K * timeouty))
Parameters:
Estimate Std. Error t value Pr(>|t|)
Linf 88.048090 3.519198 25.02 1.20e-13
K 0.084548 0.008202 10.31 3.35e-08
Residual standard error: 1.388 on 15 degrees of freedom
Correlation of Parameter Estimates:
Linf
K -0.91
Number of iterations to convergence: 4
Achieved convergence tolerance: 3.312e-07
which closely match the results in Box 5.6.