classWebsite <- "http://www.math.montana.edu/~jimrc/classes/stat505/"

 books <- read.table(paste(classWebsite,"data/books2.data",sep=""))
 books  ## bad column names, and bad first row

 books <- read.table(paste(classWebsite,"data/books2.data",sep=""), head=T)
 books

 plot(price ~ page, data=books, pch = unclass(bind))
 legend("bottomright",  c("hard","soft"), pch=c(1,2))

  ## a single line:

fit1 <- lm(price ~ page, books)
summary(fit1)
abline(fit1)

  ## two lines with the same slope
fit2pA <- lm(price ~ page + bind, books)
summary(fit2pA)
  ## bind shows up as bindp -- why?

abline(fit2pA, col="grey")  ## which level of bind?
abline(37.166 -20.36,   0.091,  col="grey")  ## adjusted for ______

  ## alternative fit showing the two intercepts:
fit2pB <- lm(price ~ page + bind + 0, books) ## +0 means no intercept (so does -1)
summary(fit2pB)                              ## Rsquared is now nonsense
abline( coef = coef(fit2pB)[c(3,1)], col="grey")

  ## two lines with different slopes:
fit2nonpA <- lm(price ~ page + bind + page:bind, books)
         ## shorthand:  lm(price ~ page * bind, books)
summary(fit2nonpA)
        ## 90% CI for change in slope using info from summary
          ## on 35 degrees of freedom
 round(  0.003704 +c(-1,1)* qt(.95, 35) *  0.024763, 3)  ## 
  ##( -0.038, 0.046)

abline(fit2nonpA, col=2, lty = 4)  ## for which group?

  ## alternative fit to add adjustments for you

fit2nonpB <- lm(price ~ 0 + bind + page:bind, books)
 summary(fit2nonpB)  ## can't use this Rsquared!
 
abline(coef=coef(fit2nonpB)[c(2,4)],  col=2, lty = 4)  ## for which group?