source("bioxm.R") library(multtest) library(lattice) library(scatterplot3d) sink() t <- as.matrix(getDataInReport()); t1 <- t[,(2:3)] X <- t(sapply(1:nrow(t1), function(r) c(as.numeric(unlist(strsplit(t1[r,], ';')))))) rownames(X) <- t[,1] num.group1 <- length(unlist(strsplit(t1[1,1],';'))) num.group2 <- length(unlist(strsplit(t1[1,2],';'))) cl <- c(rep(0, num.group1), rep(1, num.group2)) means.group1 = apply(X[,1:num.group1], 1, mean) means.group2 = apply(X[,(num.group1+1):(num.group1+num.group2)], 1, mean) means.group1 = round(means.group1, digits=2) means.group2 = round(means.group2, digits=2) resT <- mt.maxT(X, classlabel=cl, B=200) ord <- order(resT$index) adjp.BH <- mt.rawp2adjp(rawp=resT$rawp[ord], proc=c("BH")) stat <- round(resT$teststat, digits=2) res <- cbind(means.group1[ord], means.group2[ord], stat, resT$rawp, adjp.BH$adjp[ord, "BH"]) colnames(res) <- c("means.group1", "means.group2", "statistic", "rawp", "adjp.BH") res <- res[names(sort(res[,'adjp.BH'], decreasing=FALSE)), ] write.table(res, sep="\t", quote=FALSE, file=outfile) labels <- rownames(res[1:3,]) labels <- strsplit(labels, ":") xlab <- labels[[1]][3] ylab <- labels[[2]][3] zlab <- labels[[3]][3] scatterplot3d(t(X[rownames(res[1:3,]),]), xlab=xlab, ylab=ylab, zlab=zlab, color=c(rep("darkgreen", num.group1), rep("red", num.group2)), font.lab=2, cex.lab=0.6, cex.axis = 0.6, cex=0.8, pch=21, type='h' ) #color<-c(rep("green", num.group1), rep("red", num.group2))