pottery <- read.table("pottery.dat",header=TRUE)
 pottery[1:3,]
 x <- pottery[,1:5,]
 g <- pottery[,6]
 datos <- split(x,g)
 res <- lapply(datos,cov.wt)
 res

 names(res)
 res[[1]]$center
 hotelling <- function(res) {
  pobs  <- names(res)
  npobs <- length(pobs)
  for (i in 1:(npobs-1)) {
      for (j in (i+1):npobs) {
         xbar <- res[[i]]$center ; ybar <- res[[j]]$center
         nx   <- res[[i]]$n.obs  ; ny   <- res[[j]]$n.obs
         S    <- (nx*res[[i]]$cov + ny*res[[j]]$cov) / (nx+ny)
         invS <- solve(S)
         glnum <- d <- length(xbar)
         glden <- nx + ny - d - 1
         constante <- ((nx*ny)*glden) / (d*(nx + ny)^2)
         f <- constante*t(xbar - ybar)%*%invS%*%(xbar-ybar)
         cat("\nContraste igualdad medias ",pobs[i]," y ",pobs[j],"\n")
         cat("F = ",f," gl = ",glnum,glden," p-value = ",1-pf(f,glnum,glden))
         }
    }
}

 hotelling(res)