download code
## Chunk 1
library("ALL")
data("ALL")
## Chunk 2
bcell = grep("^B", as.character(ALL$BT))
moltyp = which(as.character(ALL$mol.biol) 
    %in% c("NEG", "BCR/ABL"))
ALL_bcrneg = ALL[, intersect(bcell, moltyp)]
ALL_bcrneg$mol.biol = factor(ALL_bcrneg$mol.biol)
## Chunk 3
library("genefilter")
ALLfilt_bcrneg = nsFilter(ALL_bcrneg, var.cutoff=0.5)$eset
  (solution chunk)
## Chunk 5
library("GSEABase")
gsc = GeneSetCollection(ALLfilt_bcrneg, 
                  setType=KEGGCollection())
Am = incidence(gsc)
dim(Am)
## Chunk 6
nsF = ALLfilt_bcrneg[colnames(Am),]
## Chunk 7
rtt = rowttests(nsF, "mol.biol")
rttStat = rtt$statistic
## Chunk 8
selectedRows = (rowSums(Am)>10)
Am2 = Am[selectedRows, ]
## Chunk 9
tA = as.vector(Am2 %*% rttStat)
tAadj = tA/sqrt(rowSums(Am2))
names(tA) = names(tAadj) = rownames(Am2)
## Chunk 10
stopifnot(sum(tAadj<(-8))==1L)
## Chunk 11
qqnorm(tAadj)
## Chunk 12
library("KEGG.db")
smPW = tAadj[tAadj < (-8)]
pwName = KEGGPATHID2NAME[[names(smPW)]]
pwName
## Chunk 13
KEGGmnplot(names(smPW), nsF, "hgu95av2", nsF$"mol.biol",
  pch=16, col="darkblue")
## Chunk 14
sel = as.integer(nsF$mol.biol)
KEGG2heatmap(names(smPW), nsF, "hgu95av2", 
  col=colorRampPalette(c("white", "darkblue"))(256),
  ColSideColors=c("black", "white")[sel])
## Chunk 15
stopifnot(nlevels(nsF$mol.biol)==2L)
## Chunk 16
library(Category)
set.seed(123)
NPERM = 1000
pvals = gseattperm(nsF, nsF$mol.biol, Am2, NPERM)
pvalCut = 0.025
lowC  = names(which(pvals[, 1]<=pvalCut))
highC = names(which(pvals[, 2]<=pvalCut))
## Chunk 17
## stopifnot(identical(1L,length(lowC)))
## Chunk 18
fmt = function(x) as.vector(mapply(paste, names(x), x, sep=": "))
## Chunk 19
## getPathNames(lowC)
## Chunk 20
fmt(getPathNames(lowC))
## Chunk 21
## getPathNames(highC)
## Chunk 22
fmt(getPathNames(highC))[1:5]
  (solution chunk)
  (solution chunk)
  (solution chunk)
  (solution chunk)
## Chunk 27
EGtable = toTable(hgu95av2ENTREZID[featureNames(nsF2)])
entrezUniv = unique(EGtable$gene_id)
chrMat = MAPAmat("hgu95av2", univ=entrezUniv)
rSchr = rowSums(chrMat)
  (solution chunk)
  (solution chunk)
  (solution chunk)
## Chunk 31
Ams = Am2[union(lowC, highC),]
Amx = Ams %*% t(Ams)
minS = outer(diag(Amx), diag(Amx), pmin)
overlapIndex = Amx/minS
## Chunk 32
library("lattice")
myFun = function(x) {
  dd.row = as.dendrogram(hclust(dist(x)))
  row.ord = order.dendrogram(dd.row)
  
  dd.col = as.dendrogram(hclust(dist(t(x))))
  col.ord = order.dendrogram(dd.col)
  colnames(x) = sapply(getPathNames(colnames(x)), "[", 1L)
  levelplot(x[row.ord,col.ord],  scales = list(x = list(rot = 90)),
    xlab="", ylab="", main="overlapIndex",
    col.regions=colorRampPalette(c("white", "darkblue"))(256))
}
print(myFun(overlapIndex))
## Chunk 33
ijord = function(m, i) {
  m[lower.tri(m, diag=TRUE)] = NA
  o=order(m, decreasing=TRUE, na.last=TRUE)[i]
  cbind((o-1L)%/%nrow(m)+1L, (o-1L)%%nrow(m)+1L)
}  
nm = getPathNames(colnames(overlapIndex))
pathnames = as.vector(mapply(paste, names(nm), sapply(nm, "[", 1L), sep=": "))
k=ijord(overlapIndex, 1:2)
  (solution chunk)
## Chunk 35
P04512 = Ams["04512",]
P04510 = Ams["04510",]
lm1 = lm(rttStat ~ P04512)
summary(lm1)$coefficients
lm2 = lm(rttStat ~ P04510)
summary(lm2)$coefficients
lm3 = lm(rttStat ~ P04510+P04512)
summary(lm3)$coefficients
## Chunk 36
stopifnot(
coefficients(summary(lm1))["P04512", "Pr(>|t|)"]<0.05,
coefficients(summary(lm2))["P04510", "Pr(>|t|)"]<0.05,
coefficients(summary(lm3))["P04512", "Pr(>|t|)"]>0.05,
coefficients(summary(lm3))["P04510", "Pr(>|t|)"]<0.05)
## Chunk 37
P04512.Only = ifelse(P04512 != 0 & P04510 == 0, 1, 0)
P04510.Only = ifelse(P04512 == 0 & P04510 != 0, 1, 0)
Both        = ifelse(P04512 != 0 & P04510 != 0, 1, 0)
lm4 = lm(rttStat ~ P04510.Only + P04512.Only + Both) 
## Chunk 38
## summary(lm4)
  (solution chunk)