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 Box 3.11 Vignetter. The ANCOVA is fit with 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 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)
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