// library(ascii) // setwd("c://aaaWork//web//fishR//bookex//AIFFD//Box6_6") // Sweave("Box6_6a.Rnw",driver=RweaveAsciidoc) AIFFD Box 6.6 Vignette ======================= :Author: Derek H. Ogle :Email: dogle@northland.edu :Date: 29-June-2009 :Revision: 1 == Required Packages and Setting Options ---- > library(NCStats) # view, fit.plot > library(car) # Anova ---- == Preparing Data The data for this box was entered into an external tab-delimited text file which is read into R as follows: ---- > setwd("c://aaaWork//web//fishR//bookex//AIFFD//Box6_6") > d <- read.table("box6_6.txt",header=TRUE) ---- It is important to note the structure of this data file as all of the catch data is entered into one column/vector and the lake to which that catch corresponds is listed in a separate vector. This "looks" different then the way the data are presented at the top of Box 6.6 (though it is NOT different from the way the data were entered into SAS). The structure and a view of random rows of the data frame can be seen with, ---- > str(d) 'data.frame': 17 obs. of 3 variables: $ age : int 1 2 3 4 5 6 7 8 9 1 ... $ catch: int 433 818 243 67 48 5 30 42 22 305 ... $ lake : Factor w/ 2 levels "A","B": 1 1 1 1 1 1 1 1 1 2 ... > view(d) age catch lake 1 1 433 A 3 3 243 A 11 2 491 B 14 5 30 B 15 6 49 B 17 8 6 B ---- Also, note that the catch curve methods below will require the natural log of the catches, ---- > d$logcatch <- log(d$catch) ---- == ANCOVA to Determine if Z Differs Between Lakes The method described in Box 6.6 is the same as the ANCOVA discussed in the link:../Box3_11/Box3_11.html[Box 3.11 Vignetter]. The ANCOVA is fit with +*[red]#lm()#*+ (for "linear model") where the explanatory variable side (right-hand-side) of the linear model formula uses a notation where the two explanatory variables appear to be multiplied together. In R, this apparent multiplication notation is a short-hand method of telling R to include each main effect and the interaction between the two explanatory variables. In the example below, the +*age*lake*+ portion of the formula is equivalent to saying +*age*+ + +*lake*+ + +*age:lake*+. Also note the use of the +*[red]#subset=#*+ argument to restrict the model to only age-2 and older fish (the authors in Box 6.6 note that fish age-2 and older are fully recruited to the gear). Thus, the ANOVA model using this short-hand notation is fit with, ---- > lm1 <- lm(logcatch ~ age * lake, data = d, subset = age >= 2) ---- The type-I SS are found with, ---- > anova(lm1) Analysis of Variance Table Response: logcatch Df Sum Sq Mean Sq F value Pr(>F) age 1 20.0825 20.0825 26.0612 0.0003414 lake 1 0.3946 0.3946 0.5120 0.4891677 age:lake 1 0.6598 0.6598 0.8563 0.3746383 Residuals 11 8.4765 0.7706 ---- the type-III SS are found with, ---- > Anova(lm1, type = "III") Anova Table (Type III tests) Response: logcatch Sum Sq Df F value Pr(>F) (Intercept) 101.361 1 131.5364 1.85e-07 age 21.113 1 27.3984 0.0002792 lake 0.263 1 0.3412 0.5708985 age:lake 0.660 1 0.8563 0.3746383 Residuals 8.477 11 ---- and the type-II SS are found with, ---- > Anova(lm1, type = "II") Anova Table (Type II tests) Response: logcatch Sum Sq Df F value Pr(>F) age 20.4650 1 26.5575 0.0003166 lake 0.3946 1 0.5120 0.4891677 age:lake 0.6598 1 0.8563 0.3746383 Residuals 8.4765 11 ---- Because the interaction term is the last term in this model all three types of SS provide the same p-value leading to the result that the interaction is insignificant. An insignificant interaction indicates that the slope does not differ between the groups. As the slope is directly related to the estimated instantaneous mortality rate (Z) this also suggests that there is no difference in slopes between the two lakes. A plot of the fitted models, ---- > fit.plot(lm1) ---- image::Box6_6a-010.png[] The overall estimate of the common mortality rate can be found by fitting a model to the combined data set (i.e., ignoring +*lake*+), ---- > lm2 <- lm(logcatch ~ age, data = d, subset = age >= 2) > anova(lm2) Analysis of Variance Table Response: logcatch Df Sum Sq Mean Sq F value Pr(>F) age 1 20.0825 20.0825 27.392 0.0001613 Residuals 13 9.5309 0.7331 > coef(lm2) (Intercept) age 6.7901714 -0.5320886 > confint(lm2) 2.5 % 97.5 % (Intercept) 5.538711 8.0416322 age -0.751722 -0.3124552 ----