// library(ascii)
// setwd("c://aaaWork//web//fishR//bookex//AIFFD//Box5_4")
// Sweave("Box5_4a.Rnw",driver=RweaveAsciidoc)








AIFFD Box 5.4  Vignette
=======================
:Author:    Derek H. Ogle
:Email:     dogle@northland.edu
:Date:      17-June-2009
:Revision:  2

.Author Comment
****
I have provided a more thorough analysis of how to fit von Bertalanffy models link:../../../gnrlex/VonBertalanffy/VonBertalanffy.pdf[here].  You may want to read that document before continuing with the description below for Box 5.4.
****

== Required Packages
----
> library(FSA)     # growmodel.sim

----

== Preparing Data
The working directory is set, the data file (i.e., link:box5_4.txt[]) is read, and the structure is observed with,
----
> setwd("c://aaaWork//web//fishR//bookex//AIFFD//Box5_4")
> d <- read.table("box5_4.txt", header = TRUE)
> str(d)
'data.frame':	95 obs. of  2 variables:
 $ tl : int  388 418 438 428 539 432 444 421 438 419 ...
 $ age: int  4 4 4 5 10 4 7 4 4 4 ... 
----

== Getting Initial Values for the Model Parameters
As noted in Box 5.4, the computer algorithms for fitting non-linear models are iterative and require starting values for the parameters.  An alternative to what was described in Box 5.4 for identifying starting values for the traditional von Bertalanffy equation is to plot the data and use interactive graphics to visually "fit" the von Bertalanffy model.  The parameters from this process can then serve as the starting values for the more objective non-linear modeling algorithm.  The +*[red]#growmodel.sim()#*+ function in the +*FSA*+ package provides a mechanism for visually fitting the von Bertalanffy model.  For this purpose, +*[red]#growmodel.sim()#*+ has three arguments,

- +*[red]#model=#*+: a string that identifies which model and parameterization to use (+*[red]#model="vb"#*+ is used for the traditional von Bertalanffy model used in Box 5.4),
- +*[red]#x=#*+: a formula of the form +*[red]#len ~ age#*+ where len and age generically represent the variables containing the lengths and ages,
- +*[red]#data=#*+: a data frame that contains the length and age variables.

Alternatively, if the data frame that contains the length and age variables is attached then the vector of observed ages can be sent to the +*[red]#x=#*+ argument and the vector of lengths to the +*[red]#y=#*+ argument.

The starting values for the spotted sucker data can be obtained by manipulating the slider bars on the graphic produced with,
----
> growmodel.sim("vb", tl ~ age, data = d)

----
image::sv.png[Demonstration of Finding Von Bertalanffy Starting Values]

From this process, I (your fitting may have resulted in different values) obtained starting values of 490 for Linfinity, 0.4 for K, and 0 for t0.  For ease in the non-linear model fitting procedure described next these values should be entered into a saved list object,
----
> sv <- list(Linf = 490, K = 0.4, t0 = 0)

----

== Fitting the Von Bertalanffy Model
The von Bertalanffy model can be placed into a formula object as follows,
----
> vbmdl <- tl ~ Linf * (1 - exp(-K * (age - t0)))

----

where you must make sure that +*[red]#tl#*+ and +*[red]#age#*+ are the exact names of the variables in your data frame and +*[red]#Linf#*+, +*[red]#K#*+, and +*[red]#t0#*+ 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 +*[red]#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 von Bertalanffy model is fit to the spotted sucker data with,
----
> ssvb <- nls(vbmdl, start = sv, data = d)

----

The coefficient estimates and the correlations among coefficient estimates are extracted with,
----
> summary(ssvb, correlation = TRUE)
Formula: tl ~ Linf * (1 - exp(-K * (age - t0)))

Parameters:
     Estimate Std. Error t value Pr(>|t|)
Linf 516.2333    48.7987  10.579   <2e-16
K      0.1900     0.1194   1.591    0.115
t0    -4.5422     3.4083  -1.333    0.186

Residual standard error: 18.11 on 92 degrees of freedom

Correlation of Parameter Estimates:
   Linf  K    
K  -0.99      
t0 -0.96  0.99

Number of iterations to convergence: 31 
Achieved convergence tolerance: 9.94e-06  
----

which closely match the results in Box 5.4

WARNING: One should pay very close attention to the correlations among parameter estimates shown in the results above.  These values illustrate VERY strong correlations among the parameters which hinders interpretation as many other parameter triplets would provide nearly the same model fit (this is illustrated with the lack of change in SS for the last 15 or so iterations of the model fitting as shown in Box 5.4).

Finally, one can plot the model fit with,
----
> plot(tl ~ age, data = d, pch = 16, xlab = "Age", ylab = "Total Length (mm)", xlim = c(0, 11), ylim = c(0, 550))
> ages <- seq(0, 11, 0.1)
> preds <- predict(ssvb, data.frame(age = ages))
> lines(preds ~ ages, lwd = 2, col = "red")

----
image::Box5_4a-009.png[]

WARNING: The predicted length for age-0 fish suggests issues with the data or the use of this model with these data.

== Comparison Among Groups
The authors of Chapter 5 hint at comparisons of von Bertalanffy model parameters between groups but do not provide an example of such an analysis.  An example of this analysis in R can be found in my link:../../../gnrlex/VonBertalanffy/VonBertalanffy.pdf[Von Bertalanffy] vignette.