pdf("rpd8_3r.pdf") xy <- matrix( c( 1, 0.53, 2, 1.183, 3, 1.603, 4, 1.994, 5, 2.708, 6, 3.006, 7, 3.867, 8, 4.059, 9, 4.349, 10, 4.699, 11, 4.983, 12, 5.1, 13, 5.288, 14, 5.374, 1, 0.184, 2, 0.664, 3, 1.553, 4, 1.91, 5, 2.585, 6, 3.009, 7, 3.403, 8, 3.892, 9, 4.367, 10, 4.551, 11, 4.656, 12, 4.754, 13, 4.842, 14, 4.969), byrow=T, ncol=2) nobs <- nrow(xy) y <- matrix(xy[,2],ncol=1) x1 <- xy[,1] x0 <- matrix(rep(1,nobs),ncol=1) x2 <- x1*x1 x3 <- x1*x2 ybar1 <- matrix((y[1:14]+y[15:28])/2,ncol=1) ybar <- rbind(ybar1,ybar1) x <- cbind(x0,x1,x2) betahat <- solve(t(x) %*% x) %*% t(x) %*% y yhat <- x %*% betahat sspe <- t(y-ybar) %*% (y-ybar) dfpe <- nrow(y)-nrow(ybar1) mspe <- sspe/dfpe sspe dfpe mspe sslf <- t(ybar-yhat) %*% (ybar-yhat) dflf <- nrow(ybar1)-ncol(x) mslf <- sslf/dflf sslf dflf mslf flf <- mslf/mspe fcritlf <- qf(.95,dflf,dfpe) probflf <- 1-pf(flf,dflf,dfpe) print(cbind(flf,fcritlf,probflf)) vec <- seq(0,15,0.1) plot(x1,y,xlab="Days",ylab="Density of Algae",xlim=c(0,15),ylim=c(0,8)) lines(vec,betahat[1,]+betahat[2,]*vec+betahat[3,]*vec^2) dev.off()