X<-array(rnorm(10*50*3),c(10,50,3))

apply(X,1,sum)
apply(X,2,sum)
apply(X,3,sum)

apply(X,c(1,2),sum)

##############################

n <- 17; fac <- factor(rep(1:3, length = n), levels = 1:5)
table(fac)
tapply(1:n, fac, sum)
tapply(1:n, fac, sum, simplify = FALSE)
tapply(1:n, fac, range)
tapply(1:n, fac, quantile)

n<-2000
X<-rnorm(n)
Z<-sample(c(1:4),size=n,rep=T)
Z1<-sample(c(1:3),size=n,rep=T)
Z2<-sample(c(1:2),size=n,rep=T)
Y<-Z1+Z2+X

table(Z)

tapply(X,Z,summary)
tapply(X,Z,mean)
tapply(X,Z,var)
tapply(X,X1,mean)
aggregate(X,list(Z),var)

my.data<-data.frame(X,Z,Z1,Z2,Y)
aggregate(Y~Z,data=my.data,mean)

aggregate(Y~Z1+Z2,data=my.data,mean)

tapply(Y,list(Z1,Z2),mean)

###################################################################
aggregate(PEPMNK1R~HSSEX,data=NHnames,mean,na.rm=T)

aggregate(PEPMNK1R~SMOKE+HSSEX,data=NHnames,mean,na.rm=T)
aggregate(PEPMNK1R~SMOKE+HSSEX,data=NHnames,summary)

unique(NHnames$HSSEX)

table(NHnames$HSSEX,useNA="ifany")