Dear List, why does this not work? df <- data.frame(var1 = c(3,2,1), var2 = c(6,5,4), var3 = c(9,8,7), fac = c('A', 'A', 'B')) tapply(cbind(df$var1, df$var2, df$var3), df$fac, mean) Thank you, Stefan