Below are the solutions to these exercises on “Non-Linear Model in R.”

############### # # # Exercise 1 # # # ############### # load data file Mussel <- read.csv(file.choose()) mussel.nls1 <- nls(SPECIES ~ a * AREA^b, start = list(a = 0.1, b = 1), data = Mussel) summary(mussel.nls1)

## ## Formula: SPECIES ~ a * AREA^b ## ## Parameters: ## Estimate Std. Error t value Pr(>|t|) ## a 0.8584 0.2769 3.100 0.00505 ** ## b 0.3336 0.0350 9.532 1.87e-09 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.733 on 23 degrees of freedom ## ## Number of iterations to convergence: 17 ## Achieved convergence tolerance: 1.037e-06

############### # # # Exercise 2 # # # ############### plot(resid(mussel.nls1) ~ fitted(mussel.nls1))

############### # # # Exercise 3 # # # ############### SSpower <- selfStart(~ a*x^b, function(mCall, data, LHS) { xy <- sortedXyData(mCall[["x"]], LHS, data) if(nrow(xy) < 3) { stop("Too few distinct x values to fit a power function") } z <- xy[["y"]] xy[["logx"]] <- log(xy[["x"]]) xy[["logy"]] <- log(xy[["y"]]) aux <- coef(lm(logy ~ logx, xy)) pars <- c(exp(aux[[1]]), aux[[2]]) setNames(pars, mCall[c("a", "b")]) }, c("a", "b")) # generate the nls model mussel.nls2 <- nls(SPECIES ~ SSpower(AREA, a, b), data = Mussel) summary(mussel.nls2)

## ## Formula: SPECIES ~ SSpower(AREA, a, b) ## ## Parameters: ## Estimate Std. Error t value Pr(>|t|) ## a 0.8584 0.2769 3.100 0.00505 ** ## b 0.3336 0.0350 9.532 1.87e-09 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.733 on 23 degrees of freedom ## ## Number of iterations to convergence: 5 ## Achieved convergence tolerance: 4.427e-06

############### # # # Exercise 4 # # # ############### mussel.nls2 <- nls(SPECIES ~ SSasymp(AREA, a, b, c), data = Mussel) summary(mussel.nls2)

## ## Formula: SPECIES ~ SSasymp(AREA, a, b, c) ## ## Parameters: ## Estimate Std. Error t value Pr(>|t|) ## a 24.4114 1.6644 14.667 7.71e-13 *** ## b 4.9563 1.4244 3.479 0.00213 ** ## c -8.8138 0.2482 -35.512 < 2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.719 on 22 degrees of freedom ## ## Number of iterations to convergence: 0 ## Achieved convergence tolerance: 7.122e-07

# so which one fits the data better? # note we cannot use anova() for comparison because the models are specified differently ############### # # # Exercise 5 # # # ############### AIC(mussel.nls1, mussel.nls2)

## df AIC ## mussel.nls1 3 125.1312 ## mussel.nls2 4 125.7644

# The power function appears to be the best ############### # # # Exercise 6 # # # ############### # plot in one graph plot(SPECIES ~ AREA, data = Mussel, pch = 16, xlab = (expression(paste("Mussel bed area ", "(", m^2, ")"))), ylab = "Species Richness") ############### # # # Exercise 7 # # # ############### x <- seq(min(Mussel$AREA), max(Mussel$AREA), l=1000) # 1000 steps points(x, predict(mussel.nls1, data.frame(AREA = x)), type = "l",lty = 2) # large dashed line (power function) points(x, predict(mussel.nls2, data.frame(AREA = x)), type = "l",lty = 3) # short dashed line (Asymptotic function)

## Leave a Reply