Y <- matrix(c(87,86,94,96,84,94,97,93,96,97,100,98),byrow=T,ncol=1) XF <- matrix(c( 1, 89, 0, 0, 0, 0, 1, 82, 0, 0, 0, 0, 1, 88, 0, 0, 0, 0, 1, 94, 0, 0, 0, 0, 0, 0, 1, 89, 0, 0, 0, 0, 1, 90, 0, 0, 0, 0, 1, 91, 0, 0, 0, 0, 1, 92, 0, 0, 0, 0, 0, 0, 1, 89, 0, 0, 0, 0, 1, 99, 0, 0, 0, 0, 1, 84, 0, 0, 0, 0, 1, 87), byrow=T,ncol=6) XM <- matrix(c( 1, 0, 0, 89, 1, 0, 0, 82, 1, 0, 0, 88, 1, 0, 0, 94, 0, 1, 0, 89, 0, 1, 0, 90, 0, 1, 0, 91, 0, 1, 0, 92, 0, 0, 1, 89, 0, 0, 1, 99, 0, 0, 1, 84, 0, 0, 1, 87), byrow=T, ncol=4) XR <- matrix(c( 1, 89, 1, 82, 1, 88, 1, 94, 1, 89, 1, 90, 1, 91, 1, 92, 1, 89, 1, 99, 1, 84, 1, 87), byrow=T, ncol=2) I_n <- diag(nrow(Y)) PF <- XF %*% solve(t(XF) %*% XF) %*% t(XF) PM <- XM %*% solve(t(XM) %*% XM) %*% t(XM) PR <- XR %*% solve(t(XR) %*% XR) %*% t(XR) SSE_F <- t(Y) %*% (I_n - PF) %*% Y; dfE_F <- nrow(Y) - ncol(XF) SSE_M <- t(Y) %*% (I_n - PM) %*% Y; dfE_M <- nrow(Y) - ncol(XM) SSE_R <- t(Y) %*% (I_n - PR) %*% Y; dfE_R <- nrow(Y) - ncol(XR) cat("Full Model SSE_F =", SSE_F, "dfE_F =", dfE_F, "\n") cat("Middle Model SSE_M =", SSE_M, "dfE_M =", dfE_M, "\n") cat("Reduced Model SSE_R =", SSE_R, "dfE_R =", dfE_R, "\n") #### Test for Equality of Slopes Full Model vs Middle Model F1 <- ((SSE_M-SSE_F)/(dfE_M-dfE_F))/(SSE_F/dfE_F) F1_05 <- qf(.95,dfE_M-dfE_F,dfE_F) probF1 <- 1-pf(F1,dfE_M-dfE_F,dfE_F) #### Test for Equality of Intercepts, Given Equal Slopes - Mid vs Reduced F2 <- ((SSE_R-SSE_M)/(dfE_R-dfE_M))/(SSE_M/dfE_M) F2_05 <- qf(.95,dfE_R-dfE_M,dfE_M) probF2 <- 1-pf(F2,dfE_R-dfE_M,dfE_M) ### Joint Test for Equality of Slopes and Intercepts - Full vs Reduced F3 <- ((SSE_R-SSE_F)/(dfE_R-dfE_F))/(SSE_F/dfE_F) F3_05 <- qf(.95,dfE_R-dfE_F,dfE_F) probF3 <- 1-pf(F3,dfE_R-dfE_F,dfE_F) cat("Full vs Middle: F1 =", F1, "F1_05 =", F1_05, "probF1 =", probF1, "\n") cat("Middle vs Reduced: F2 =", F2, "F2_05 =", F2_05, "probF2 =", probF2, "\n") cat("Full vs Reduced: F3 =", F3, "F3_05 =", F3_05, "probF3 =", probF3, "\n") ##### Bartlett's Test for equal variances e_F <- (I_n - PF) %*% Y SSE1 <- t(e_F[1:4,1]) %*% e_F[1:4,1]; dfE1 <- 4-2; s2_1 <- SSE1/dfE1; SSE2 <- t(e_F[5:8,1]) %*% e_F[5:8,1]; dfE2 <- 4-2; s2_2 <- SSE2/dfE2; SSE3 <- t(e_F[9:12,1]) %*% e_F[9:12,1]; dfE3 <- 4-2; s3_2 <- SSE3/dfE3; dfE <- dfE1+dfE2+dfE3 MSE <- (SSE1 + SSE2 + SSE3)/dfE t <- 3 C <- 1 + (1/(3*(t-1)))*(1/dfE1 + 1/dfE2 + 1/dfE3 - 1/dfE) B <- (1/C)*(dfE*log(MSE) - (dfE1*log(s2_1)+dfE2*log(s2_2)+dfE3*log(s3_2))) X2B_05 <- qchisq(.95,t-1); probB <- 1-pchisq(B,t-1); cat("Bartlett's Test: B=", B, "X2_05 =", X2B_05, "probB =", probB, "\n")