###################################################
### chunk number 1: initialize
###################################################
library(safe)
library(multtest)
library(hu6800.db)


###################################################
### chunk number 2: 
###################################################
data(golub)
dimnames(golub)[[1]] <- golub.gnames[,3]


###################################################
### chunk number 3: 
###################################################
table(golub.cl)


###################################################
### chunk number 4: 
###################################################
C.mat <- getCmatrix(gene.list = as.list(hu6800PATH), as.matrix = TRUE,
                    present.genes = golub.gnames[,3], min.size = 10)
dimnames(C.mat)[[2]] <- paste("KEGG:",dimnames(C.mat)[[2]],sep="") 


###################################################
### chunk number 5: 
###################################################
set.seed(12345)
results <- safe(golub, golub.cl, platform = "hu6800",annotate = "KEGG", min.size = 10)


###################################################
### chunk number 6: 
###################################################
results


###################################################
### chunk number 7: 
###################################################
gene.results(results,cat.name="KEGG:00860")


###################################################
### chunk number 8: 
###################################################
y.vec <- c("ALL","AML")[1+golub.cl]
results2 <- safe(golub, y.vec, C.mat, Pi.mat = 1)
results3 <- safe(golub, golub.cl, C.mat, local="t.Welch", Pi.mat = 1)
round(cbind(Student1 = results@local.stat[1:3],
            Student2 = results2@local.stat[1:3],
            Welch = results3@local.stat[1:3]),3)


###################################################
### chunk number 9: 
###################################################
y.vec <- rep(1:19,2)*rep(c(-1,1),each=19)
y.vec
results2 <- safe(golub, y.vec, C.mat, local="t.paired",
                 Pi.mat = 1)


###################################################
### chunk number 10: 
###################################################
y.vec <- rexp(38)
cens <- rep(0:1,c(30,8))
results2 <- safe(golub, y.vec, C.mat, local="z.COXPH", Pi.mat = 1, 
                 args.local = list(censor=cens))


###################################################
### chunk number 11: 
###################################################
local.Wilcoxon<-function(X.mat,y.vec, ...){
  return(function(data,trt = (y.vec == 1)) {
    return(as.numeric(trt %*% apply(data,1,rank)))
  })
}
results2 <-  safe(golub, golub.cl, C.mat, Pi.mat = 1,
                  local="Wilcoxon")


###################################################
### chunk number 12: 
###################################################
cbind(Student1 = round(results@local.stat[1:3],3),
      Rank.Sum=results2@local.stat[1:3])


###################################################
### chunk number 13: 
###################################################
set.seed(12345)
results2 <- safe(golub, golub.cl, C.mat, global="Fisher",
                 args.global = list(one.sided=F,genelist.length=200))


###################################################
### chunk number 14: 
###################################################
results2


###################################################
### chunk number 15: 
###################################################
1-phyper(12-1, 70, 3051-12, 200)


###################################################
### chunk number 16: 
###################################################
set.seed(12345)
results <- safe(golub, golub.cl, C.mat, error = "FDR.YB", alpha = 0.25)


###################################################
### chunk number 17: 
###################################################
results


###################################################
### chunk number 18: 
###################################################
set.seed(12345)
results2 <- safe(golub, golub.cl, C.mat, method = "bootstrap", 
                 error = "FDR.BH")


###################################################
### chunk number 19: 
###################################################
results2


###################################################
### chunk number 20: 
###################################################
results2 <- safe(golub, golub.cl, platform="hu6800", annotate="PFAM",
                 min.size=10,method="bootstrap")


###################################################
### chunk number 21: 
###################################################
results2


###################################################
### chunk number 22: 
###################################################
GO.list <- as.list(hu6800GO2ALLPROBES)
C.mat.CC <- getCmatrix(keyword.list = GO.list, GO.ont = "CC", 
            present.genes = dimnames(golub)[[1]], min.size = 10,
            max.size=200)
results2 <- safe(golub, golub.cl, C.mat.CC, method="bootstrap")


###################################################
### chunk number 23: plot1
###################################################
safeplot(results)


###################################################
### chunk number 24: plot2
###################################################
safeplot(results,cat.name="KEGG:00010")


###################################################
### chunk number 25: plot3
###################################################
safedag(results2,filter=1)


###################################################
### chunk number 26: plot4
###################################################
safedag(results2, filter=1,top="GO:0044428")