platform="linux"
cosadir="COSADIR"
source("COSADIR/r_cosa.q")

#####
set.seed(1)
n0<-150      #  number of attributes defining the clustering
n<-10000     #  number of attributes total 
N<-100       #  number of objects
N0<-15       #  number of objects in the small cluster

Y1<-matrix( rnorm(n*(N-N0)),nrow=N-N0,ncol=n)        # simulated attributes
Y22<-matrix( rnorm( N0*(n-n0)), nrow=N0,ncol=n-n0)   # for both groups
                                                     # except the N0*n0
                                                     # attribute values 
                                                     # defining the clustering

set.seed(1)
Y21<-matrix( rnorm(N0*n0,1.5,.2), nrow=N0,ncol=n0)   # change in mean and sd
Y<-rbind(Y1,cbind(Y21,Y22))
Y<-t( t(Y)-apply(Y,2,mean))
Y<-t( t(Y)/apply(Y,2,sd))
dcosa.sdmn<-cosadist(Y)

set.seed(1)
Y21<-matrix( rnorm(N0*n0,0,.2), nrow=N0,ncol=n0)     # change in sd only 
Y<-rbind(Y1,cbind(Y21,Y22))
Y<-t( t(Y)-apply(Y,2,mean))
Y<-t( t(Y)/apply(Y,2,sd))
dcosa.sd<-cosadist(Y)
                                                                             
set.seed(1)
Y21<-matrix( rnorm(N0*n0,1.5,1), nrow=N0,ncol=n0)    # change in mean only
Y<-rbind(Y1,cbind(Y21,Y22))
Y<-t( t(Y)-apply(Y,2,mean))
Y<-t( t(Y)/apply(Y,2,sd))
dcosa.mn<-cosadist(Y)


par(mfrow=c(1,3))
groups<-c( rep(1,N-N0),rep(2,N0))
hierclust(dcosa.sdmn,method="average",labels=as.character(groups),xlab="sdmn")
hierclust(dcosa.mn,method="average",labels=as.character(groups),xlab="mn")
hierclust(dcosa.sd,method="average",labels=as.character(groups),xlab="sd")