airperm <- read.csv("http://www.stat.ufl.edu/~winner/sta4210/mydata/airperm_woven_reg.csv", header=TRUE) attach(airperm) names(airperm) plot(airperm[,2:5]) ### Scatterplot Matrix of X1,X2,X3,Y ssto <- sum((mean_ap-mean(mean_ap))^2) ### Total Sum of Squares ########## Fit all 7 possible Models ap.mod123 <- lm(mean_ap ~ warp + weft + mass) summary(ap.mod123) anova(ap.mod123) ap.mod12 <- lm(mean_ap ~ warp + weft) #summary(ap.mod12) #anova(ap.mod12) ap.mod13 <- lm(mean_ap ~ warp + mass) #summary(ap.mod13) #anova(ap.mod13) ap.mod23 <- lm(mean_ap ~ weft + mass) #summary(ap.mod23) #anova(ap.mod23) ap.mod1 <- lm(mean_ap ~ warp) #summary(ap.mod1) #anova(ap.mod1) ap.mod2 <- lm(mean_ap ~ weft) #summary(ap.mod2) #anova(ap.mod2) ap.mod3 <- lm(mean_ap ~ mass) #summary(ap.mod3) #anova(ap.mod3) ####### Obtain SSE and SSR for each model (deviance=SSE) sse.x1 <- deviance(ap.mod1); ssr.x1 <- ssto-sse.x1 sse.x2 <- deviance(ap.mod2); ssr.x2 <- ssto-sse.x2 sse.x3 <- deviance(ap.mod3); ssr.x3 <- ssto-sse.x3 sse.x1x2 <- deviance(ap.mod12); ssr.x1x2 <- ssto-sse.x1x2 sse.x1x3 <- deviance(ap.mod13); ssr.x1x3 <- ssto-sse.x1x3 sse.x2x3 <- deviance(ap.mod23); ssr.x2x3 <- ssto-sse.x2x3 sse.x1x2x3 <- deviance(ap.mod123); ssr.x1x2x3 <- ssto-sse.x1x2x3 #### Compute Sequential Sums of Squares ssr.x1x2-ssr.x1 ### SSR(X2|X1) ssr.x1x2-ssr.x2 ### SSR(X2|X1) ssr.x1x2x3-ssr.x1x2 ### SSR(X3|X1,X2) #### Test H0: B2=B3=0 anova(ap.mod1,ap.mod123) ### Compue Coefficients of Partial Determination (ssr.x1x2-ssr.x1)/sse.x1 ### R2(Y2|1) (ssr.x1x2x3-ssr.x1x2)/sse.x1x2 ### R2(Y3|12) #### Correlation Transformation and Standardized Regression Coefficients y.corr <- (mean_ap-mean(mean_ap))/(sqrt(29)*sd(mean_ap)) x1.corr <- (warp-mean(warp))/(sqrt(29)*sd(warp)) x2.corr <- (weft-mean(weft))/(sqrt(29)*sd(weft)) x3.corr <- (mass-mean(mass))/(sqrt(29)*sd(mass)) xstar <- matrix(cbind(x1.corr,x2.corr,x3.corr),ncol=3) ystar <- matrix(y.corr) bstar <- solve(t(xstar) %*% xstar) %*% t(xstar) %*% ystar bstar ##### Regression coefficient estimates for all 2 and 3 variable models summary(ap.mod12) summary(ap.mod13) summary(ap.mod23) summary(ap.mod123)