"d.outlier" <- function(da){
# Program to detect outliers
nr = dim(da)[1]
nc = dim(da)[2]
cm=matrix(colMeans(da),1,nc)
o=matrix(rep(1,nr),nr,1)
dev=da-kronecker(o,cm)
dev=as.matrix(dev)      
s=t(dev)%*%dev/(nr-1)             
si=solve(s)
d2=diag(dev%*%si%*%t(dev))
sd=sqrt(diag(s))
z = dev
for (i in 1:nc) {
  z[,i]= z[,i]/sd[i]
}
pairs(z)
result=cbind(z,d2)
signif(result,2)
}