# 
# Written by:
# -- 
# John L. Weatherwax                2009-04-21
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
#-----

library(alr3)
data(snake)
attach(snake)

# 2.7.3
#

m <- lm( Y ~ X - 1 ) 

summary(m)

postscript("../../WriteUp/Graphics/Chapter2/prob_7_snake.eps", onefile=FALSE, horizontal=FALSE)
plot( X, Y, xlab="water content", ylab="water yield" ) 
abline( m ) 
dev.off()

# alternativly we can compute our estimate "from scratch"
#

beta_hat_1 = sum( X * Y ) / sum( X * X )
RSS        = sum( m$residuals^2 )
n          = length( X ) 
sigma_hat  = sqrt( RSS / (n-1) )
se_beta_hat_1 = sqrt( sigma_hat^2 / sum( X*X ) )

# a 95% confidence interval for beta_1
alpha = 0.05
quant = qt( 1 - alpha/2, n-1 ) 

beta_hat_1 - quant * se_beta_hat_1
beta_hat_1 + quant * se_beta_hat_1

# 2.7.4: plot the residuals as a function of fit value y_hat:
#
y_hat <- predict(m)

postscript("../../WriteUp/Graphics/Chapter2/prob_7_res_vs_fit.eps", onefile=FALSE, horizontal=FALSE)
plot( y_hat, m$residuals, xlab="Fitted values", ylab="Residuals" )
dev.off()

sum( m$residuals )