Plan <- read.table("https://husson.github.io/livre_exo/chap13_ex20.csv", header = TRUE, sep=";",row.names=1, encoding = "latin1") for (i in 1:3) Plan[,i] <- as.factor(Plan[,i]) X <- model.matrix(~ F1+F2+F3+F1:F2+F1:F3+F2:F3, Plan, contrasts = "contr.sum") t(X)%*%X solve(t(X)%*%X) X <- model.matrix(~ F1+F2+F3, Plan, contrasts = "contr.sum") t(X)%*%X solve(t(X)%*%X) library(FactoMineR) AovSum(Y~ F1+F2+F3, Plan)