###################################################
### chunk number 1: options
###################################################
options(continue = "  ")


###################################################
### chunk number 2: getGolubData
###################################################
library(golubEsets)
library(vsn)
data(Golub_Merge)
Golub<-Golub_Merge[,1:34]
exprs(Golub)<-exprs(vsn2(Golub_Merge[,1:34]))


###################################################
### chunk number 3: getNKIData
###################################################
library(globaltest)
data(vandeVijver)


###################################################
### chunk number 4: loadKEGGandGO
###################################################
library(hu6800.db)
kegg <- as.list(hu6800PATH2PROBE)
go <- as.list(hu6800GO2ALLPROBES)


###################################################
### chunk number 5: cleanGO
###################################################
go <- lapply(go, function(x) x[!is.na(names(x)) & (names(x) != "IEA")])


###################################################
### chunk number 6: makeCellcycle
###################################################
GO.cellcycle <- go[["GO:0007049"]]
KEGG.cellcycle <- kegg[["04110"]]


###################################################
### chunk number 7: all
###################################################
gt.all <- globaltest(Golub, "ALL.AML")


###################################################
### chunk number 8: all.display
###################################################
gt.all


###################################################
### chunk number 9: cellcycle
###################################################
globaltest(Golub, "ALL.AML", GO.cellcycle)


###################################################
### chunk number 10: twocellcycle
###################################################
cellcycle <- list(go = GO.cellcycle, kegg = KEGG.cellcycle)
globaltest(Golub, "ALL.AML", cellcycle)


###################################################
### chunk number 11: kegg
###################################################
gt.kegg <- globaltest(Golub, "ALL.AML", kegg)


###################################################
### chunk number 12: kegg.display eval=FALSE
###################################################
## gt.kegg["04110"]
## gt.kegg[1:10]


###################################################
### chunk number 13: sort.by.p
###################################################
sorted <- sort(gt.kegg)
top5 <- sorted[1:5]
top5


###################################################
### chunk number 14: kegg.numbers.to.names
###################################################
library(KEGG.db)
names(top5) <- as.list(KEGGPATHID2NAME)[names(top5)]
top5


###################################################
### chunk number 15: useful
###################################################
p.value(top5)
names(top5)
result(top5)
length(top5)


###################################################
### chunk number 16: multtest
###################################################
sorted <- gt.multtest(sorted, "FWER")
sorted[1:5]


###################################################
### chunk number 17: golub_Source
###################################################
globaltest(Golub, "Source")


###################################################
### chunk number 18: golub_Source_2levels
###################################################
globaltest(Golub, "Source", levels = c("CALGB", "CCG"))


###################################################
### chunk number 19: golub_Source_1level
###################################################
globaltest(Golub, "Source", levels = "St-Jude")


###################################################
### chunk number 20: golub_Source_factor eval=FALSE
###################################################
## globaltest(Golub, "factor(Source)")


###################################################
### chunk number 21: golub_pctBlasts
###################################################
calgb <- pData(Golub)["Source"] == "CALGB"
globaltest(Golub[,calgb], "pctBlasts")


###################################################
### chunk number 22: vandeVijver_survival
###################################################
library(survival)
globaltest(vandeVijver, "Surv(TIMEsurvival, EVENTdeath)")


###################################################
### chunk number 23: vandeVijver_survival_alternatives eval=FALSE
###################################################
## globaltest(vandeVijver, "Surv(TIMEsurvival, EVENTdeath == 1)")
## globaltest(vandeVijver, "TIMEsurvival", d = "EVENTdeath")
## globaltest(vandeVijver, "TIMEsurvival", d = "EVENTdeath", event = 1)


###################################################
### chunk number 24: method
###################################################
globaltest(Golub, "ALL.AML", GO.cellcycle, method = "p", nperm = 1000)


###################################################
### chunk number 25: permutations
###################################################
gt <- globaltest(Golub, "ALL.AML", GO.cellcycle)
permutations(gt, nperm = 2000)


###################################################
### chunk number 26: golub_adjust_BM.PB
###################################################
globaltest(Golub, ALL.AML ~ BM.PB, GO.cellcycle)


###################################################
### chunk number 27: vijver_adjust_all
###################################################
globaltest(vandeVijver, Surv(TIMEsurvival, EVENTdeath) ~ NIH + ESR1 + Posnodes)


###################################################
### chunk number 28: vijver_esr1
###################################################
globaltest(vandeVijver, "ESR1")


###################################################
### chunk number 29: golub_adjust_AML.ALL
###################################################
gt <- globaltest(vandeVijver, Surv(TIMEsurvival, EVENTdeath) ~ ESR1)
fit(gt)


###################################################
### chunk number 30: makeGOstructure
###################################################
bp <- makeGOstructure(Golub, "hu6800", unreliable = "IEA")
bp


###################################################
### chunk number 31: getFocus
###################################################
focusBP <- getFocus(bp, maxatoms = 7)


###################################################
### chunk number 32: gtGO
###################################################
go60 <- gtGO(Golub, "ALL.AML", focus = focusBP, GO = bp, stopafter = 60)


###################################################
### chunk number 33: visualizeGO
###################################################
library(GOstats)
library(Rgraphviz)
sigGO <- GOGraph(names(go60), GOBPPARENTS)
sigGO <- removeNode("all", sigGO)
sigGO <- as(t(as(sigGO, "matrix")),"graphNEL")
nodes <- buildNodeList(sigGO)
focusnode <- sapply(nodes, name) %in% focusBP
names(focusnode) <- names(nodes)
nodefill <- ifelse(focusnode, "#BBBBBB", "white")
nAttrs <- list()
nAttrs$fillcolor <- nodefill
nAttrs$label <- 1:length(names(nodes))
names(nAttrs$label) <- names(nodes)
pg <- plot(sigGO, nodeAttrs = nAttrs)
x <- getNodeXY(pg)$x
y <- getNodeXY(pg)$y
ordering <- sort.list(order(-y,x))
nAttrs$label <- ordering
names(nAttrs$label) <- names(nodes)
plot(sigGO, nodeAttrs = nAttrs)
Terms <- sapply(mget(names(nodes)[sort.list(ordering)], GOTERM), Term)
names(Terms) <- NULL
legend <- data.frame(Terms)


###################################################
### chunk number 34: graph-plot
###################################################
plot(sigGO, nodeAttrs = nAttrs)


###################################################
### chunk number 35: graphlegend
###################################################
legend


###################################################
### chunk number 36: interactivity eval=FALSE
###################################################
## repeat {
##   p <- locator(n=1)
##   if (is.null(p)) break()
##   pg <- plot(sigGO, nodeAttrs = nAttrs)
##   x <- getNodeXY(pg)$x
##   y <- getNodeXY(pg)$y
##   distance <- abs(p$x-x)+abs(p$y-y)
##   idx <- which.min(distance)
##   legend("topleft",legend=c(nAttrs$label[idx],names(focusnode)[idx],
##     Term(mget(names(focusnode)[idx], GOTERM)[[1]]),bg="white"))
##   }


###################################################
### chunk number 37: sampling
###################################################
gt <- globaltest(Golub, "ALL.AML", KEGG.cellcycle)
sampled.gt <- sampling(gt)
sampled.gt


###################################################
### chunk number 38: hist eval=FALSE
###################################################
## gt <- globaltest(Golub, "ALL.AML", KEGG.cellcycle, method = "p")
## hist(gt)


###################################################
### chunk number 39: hist-plot
###################################################
gt <- globaltest(Golub, "ALL.AML", KEGG.cellcycle, method = "p")
hist(gt)


###################################################
### chunk number 40: geneplot
###################################################
gt <- globaltest(Golub, "ALL.AML", kegg)
geneplot(gt, "00561")


###################################################
### chunk number 41: geneplot_store
###################################################
myplot <- geneplot(gt, "00561", plot = FALSE)


###################################################
### chunk number 42: geneplot_changenames
###################################################
names(myplot) <- as.list(hu6800SYMBOL)[names(myplot)]
plot(myplot)


###################################################
### chunk number 43: geneplot_changenames2
###################################################
plot(myplot)


###################################################
### chunk number 44: geneplot_sort
###################################################
mysorted <- sort(myplot)
top5 <- mysorted[1:5]


###################################################
### chunk number 45: geneplot_utilities
###################################################
z.score(top5)
names(top5)
result(top5)
length(top5)


###################################################
### chunk number 46: sampleplot
###################################################
gt <- globaltest(Golub, "ALL.AML", KEGG.cellcycle)
sampleplot(gt)


###################################################
### chunk number 47: sampleplot-plot
###################################################
sampleplot(gt)


###################################################
### chunk number 48: checkerboard
###################################################
gt <- globaltest(Golub, "ALL.AML", kegg)
checkerboard(gt.kegg, "04110")


###################################################
### chunk number 49: checkerboard-plot
###################################################
checkerboard(gt, "04110")


###################################################
### chunk number 50: regressionplot
###################################################
gt <- globaltest(Golub, "ALL.AML", kegg)
regressionplot(gt, "04110", sampleid = "39")


###################################################
### chunk number 51: regressionplot2
###################################################
regressionplot(gt, "04110", sampleid = c("39","40"))


###################################################
### chunk number 52: regressionplot-plot
###################################################
regressionplot(gt, "04110", "39")


###################################################
### chunk number 53: 
###################################################
toLatex(sessionInfo())