################################################### ### chunk number 1: ################################################### #line 47 "vignettes/XDE/inst/doc/XDE.Rnw" options(width=60) ################################################### ### chunk number 2: data ################################################### #line 116 "vignettes/XDE/inst/doc/XDE.Rnw" library(XDE) data(expressionSetList) xlist <- expressionSetList class(xlist) ################################################### ### chunk number 3: validObject ################################################### #line 133 "vignettes/XDE/inst/doc/XDE.Rnw" validObject(xlist) ################################################### ### chunk number 4: covariate ################################################### #line 146 "vignettes/XDE/inst/doc/XDE.Rnw" all(sapply(xlist, function(x, label){ label %in% varLabels(x)}, label="adenoVsquamous")) ################################################### ### chunk number 5: initializeList ################################################### #line 185 "vignettes/XDE/inst/doc/XDE.Rnw" params <- new("XdeParameter", esetList=xlist, phenotypeLabel="adenoVsquamous") params ################################################### ### chunk number 6: initializeList2 eval=FALSE ################################################### ## #line 196 "vignettes/XDE/inst/doc/XDE.Rnw" ## params <- new("XdeParameter", esetList=xlist, ## phenotypeLabel="adenoVsquamous", one.delta=FALSE) ################################################### ### chunk number 7: empiricalStart eval=FALSE ################################################### ## #line 231 "vignettes/XDE/inst/doc/XDE.Rnw" ## params <- new("XdeParameter", esetList=xlist, phenotypeLabel="adenoVsquamous", one.delta=FALSE) ## empirical <- empiricalStart(xlist, phenotypeLabel="adenoVsquamous") ## firstMcmc(params) <- empirical ################################################### ### chunk number 8: burnin eval=FALSE ################################################### ## #line 239 "vignettes/XDE/inst/doc/XDE.Rnw" ## iterations(params) <- 3 ## burnin(params) <- TRUE ## fit <- xde(params, xlist) ################################################### ### chunk number 9: moreIterations eval=FALSE ################################################### ## #line 255 "vignettes/XDE/inst/doc/XDE.Rnw" ## iterations(params) <- 2 ## fit2 <- xde(params, xlist, fit) ################################################### ### chunk number 10: firstMcmc eval=FALSE ################################################### ## #line 266 "vignettes/XDE/inst/doc/XDE.Rnw" ## firstMcmc(params) <- lastMcmc(fit) ## seed(params) <- seed(fit) ## fit2 <- xde(params, xlist) ################################################### ### chunk number 11: runMoreIterations eval=FALSE ################################################### ## #line 281 "vignettes/XDE/inst/doc/XDE.Rnw" ## burnin(params) <- FALSE ## iterations(params) <- 2000 ## output(params)[c("potential", "acceptance", ## "diffExpressed", ## "nu", ##"DDelta", ## ##"delta", ## "probDelta", ## ##"sigma2", ## "phi")] <- 0 ## thin(params) <- 2 ## directory(params) <- "logFiles" ## xmcmc <- xde(params, xlist) ################################################### ### chunk number 12: savexmcmc eval=FALSE ################################################### ## #line 296 "vignettes/XDE/inst/doc/XDE.Rnw" ## ##this function should work ## postAvg <- calculatePosteriorAvg(xmcmc, NCONC=2, NDIFF=1) ## save(postAvg, file="~/madman/Rpacks/XDE/inst/logFiles/postAvg.rda") ## ##put browser in .standardizedDelta to fix tau2R, ... ## BES <- calculateBayesianEffectSize(xmcmc) ## save(BES, file="~/madman/Rpacks/XDE/inst/logFiles/BES.rda") ## save(xmcmc, file="~/madman/Rpacks/XDE/data/xmcmc.RData") ## q("no") ################################################### ### chunk number 13: loadObject ################################################### #line 311 "vignettes/XDE/inst/doc/XDE.Rnw" data(xmcmc, package="XDE") ################################################### ### chunk number 14: savedParams ################################################### #line 324 "vignettes/XDE/inst/doc/XDE.Rnw" output(xmcmc)[2:22][output(xmcmc)[2:22] == 1] ################################################### ### chunk number 15: extractc2 ################################################### #line 336 "vignettes/XDE/inst/doc/XDE.Rnw" pathToLogFiles <- system.file("logFiles", package="XDE") xmcmc@directory <- pathToLogFiles c2 <- xmcmc$c2 ################################################### ### chunk number 16: c2 ################################################### #line 342 "vignettes/XDE/inst/doc/XDE.Rnw" par(las=1) plot.ts(c2, ylab="c2", xlab="iterations", plot.type="single") ################################################### ### chunk number 17: retrieveLogs ################################################### #line 351 "vignettes/XDE/inst/doc/XDE.Rnw" getLogs <- function(object){ params <- output(object)[output(object) == 1] params <- params[!(names(params) %in% c("nu", "phi", "DDelta", "delta", "sigma2", "diffExpressed"))] names(params) } param.names <- getLogs(xmcmc) params <- lapply(lapply(as.list(param.names), function(name, object) eval(substitute(object$NAME_ARG, list(NAME_ARG=name))), object=xmcmc), as.ts) names(params) <- param.names tracefxn <- function(x, name) plot(x, plot.type="single", col=1:ncol(x), ylab=name) mapply(tracefxn, params, name=names(params)) ################################################### ### chunk number 18: bayesianEffectSize eval=FALSE ################################################### ## #line 371 "vignettes/XDE/inst/doc/XDE.Rnw" ## bayesianEffectSize(xmcmc) <- calculateBayesianEffectSize(xmcmc) ################################################### ### chunk number 19: bes ################################################### #line 379 "vignettes/XDE/inst/doc/XDE.Rnw" load(file.path(pathToLogFiles, "BES.rda")) load(file.path(pathToLogFiles, "postAvg.rda")) ################################################### ### chunk number 20: postAvgFig_conc ################################################### #line 392 "vignettes/XDE/inst/doc/XDE.Rnw" par(las=1) hist(postAvg[, "concordant"], breaks=50, xlim=c(0, 1), main="", xlab="posterior probability of concordant differential expression") ################################################### ### chunk number 21: postAvgFig_disc ################################################### #line 398 "vignettes/XDE/inst/doc/XDE.Rnw" par(las=1) hist(postAvg[, "discordant"], breaks=50, xlim=c(0, 1), main="", xlab="posterior probability of discordant differential expression") ################################################### ### chunk number 22: setPar ################################################### #line 411 "vignettes/XDE/inst/doc/XDE.Rnw" op.conc <- symbolsInteresting(rankingStatistic=postAvg[, "concordant"]) op.disc <- symbolsInteresting(rankingStatistic=postAvg[, "discordant"]) ################################################### ### chunk number 23: conc ################################################### #line 416 "vignettes/XDE/inst/doc/XDE.Rnw" par(las=1) graphics:::pairs(BES[op.conc$order, ], pch=op.conc$pch, col=op.conc$col, bg=op.conc$bg, upper.panel=NULL, cex=op.conc$cex) ################################################### ### chunk number 24: disc ################################################### #line 424 "vignettes/XDE/inst/doc/XDE.Rnw" graphics:::pairs(BES[op.disc$order, ], pch=op.disc$pch, col=op.disc$col, bg=op.disc$bg, upper.panel=NULL, cex=op.disc$cex) ################################################### ### chunk number 25: ssAlternatives eval=FALSE ################################################### ## #line 437 "vignettes/XDE/inst/doc/XDE.Rnw" ## ##t <- ssStatistic(statistic="t", phenotypeLabel="adenoVsquamous", esetList=xlist) ## tt <- rowttests(xlist, "adenoVsquamous", tstatOnly=TRUE) ## if(require(siggenes)){ ## sam <- ssStatistic(statistic="sam", phenotypeLabel="adenoVsquamous", esetList=xlist) ## } ## if(require(GeneMeta)){ ## z <- ssStatistic(statistic="z", phenotypeLabel="adenoVsquamous", esetList=xlist) ## } ################################################### ### chunk number 26: tstatConc eval=FALSE ################################################### ## #line 448 "vignettes/XDE/inst/doc/XDE.Rnw" ## graphics:::pairs(tt[op.conc$order, ], pch=op.conc$pch, col=op.conc$col, ## bg=op.conc$bg, upper.panel=NULL, cex=op.conc$cex) ## ################################################### ### chunk number 27: zConc eval=FALSE ################################################### ## #line 454 "vignettes/XDE/inst/doc/XDE.Rnw" ## graphics:::pairs(tt[op.conc$order, ], pch=op.conc$pch, col=op.conc$col, ## bg=op.conc$bg, upper.panel=NULL, cex=op.conc$cex) ## ################################################### ### chunk number 28: tstatDisc eval=FALSE ################################################### ## #line 460 "vignettes/XDE/inst/doc/XDE.Rnw" ## graphics:::pairs(tt[op.disc$order, ], pch=op.disc$pch, col=op.disc$col, ## bg=op.disc$bg, upper.panel=NULL, cex=op.disc$cex) ################################################### ### chunk number 29: pca eval=FALSE ################################################### ## #line 489 "vignettes/XDE/inst/doc/XDE.Rnw" ## tScores <- xsScores(tt, N=nSamples(xlist)) ## samScores <- xsScores(sam, nSamples(xlist)) ## zScores <- xsScores(z[, match(names(xlist), colnames(z))], N=nSamples(xlist)) ## ##Concordant differential expression, we use the combined score from the random effects model directly ## zScores[, "concordant"] <- z[, "zSco"] ################################################### ### chunk number 30: sessionInfo ################################################### #line 505 "vignettes/XDE/inst/doc/XDE.Rnw" toLatex(sessionInfo())