### R code from vignette source 'vignettes/BitSeq/inst/doc/BitSeq.Rnw' ################################################### ### code chunk number 1: BitSeq.Rnw:29-30 ################################################### options(width = 40) ################################################### ### code chunk number 2: BitSeq.Rnw:50-52 (eval = FALSE) ################################################### ## source("http://www.bioconductor.org/biocLite.R") ## biocLite("BitSeq") ################################################### ### code chunk number 3: BitSeq.Rnw:56-57 ################################################### library(BitSeq) ################################################### ### code chunk number 4: BitSeq.Rnw:86-87 ################################################### setwd(system.file("extdata",package="BitSeq")); ################################################### ### code chunk number 5: BitSeq.Rnw:98-105 ################################################### res1 <- getExpression("data-c0b0.sam", "ensSelect1.fasta", log = TRUE, MCMC_burnIn=200, MCMC_samplesN=200, MCMC_samplesSave=50, seed=47) ################################################### ### code chunk number 6: BitSeq.Rnw:112-113 (eval = FALSE) ################################################### ## hist(res1$exp$mean) ################################################### ### code chunk number 7: BitSeq.Rnw:117-118 ################################################### samples1 <- loadSamples(res1$fn) ################################################### ### code chunk number 8: BitSeq.Rnw:121-123 (eval = FALSE) ################################################### ## plot( unlist(s2["ENST00000436661",]), ## unlist(s2["ENST00000373501",])) ################################################### ### code chunk number 9: BitSeq.Rnw:134-136 ################################################### cond1Files = c(res1$fn,"data-c0b1.rpkm") cond2Files = c("data-c1b1.rpkm","data-c1b1.rpkm") ################################################### ### code chunk number 10: BitSeq.Rnw:142-145 ################################################### de1 <- getDE(list(cond1Files, cond2Files), samples=TRUE) print(de1$fn) ################################################### ### code chunk number 11: BitSeq.Rnw:149-150 ################################################### head( de1$pplr[ order(abs(0.5-de1$pplr$pplr), decreasing=TRUE ), ], 5) ################################################### ### code chunk number 12: BitSeq.Rnw:159-160 ################################################### setwd(system.file("extdata",package="BitSeq")) ################################################### ### code chunk number 13: BitSeq.Rnw:170-175 ################################################### parseAlignment( "data-c0b0.sam", outFile = "data-c0b0.prob", trSeqFile = "ensSelect1.fasta", trInfoFile = "data.tr", verbose = TRUE ) ################################################### ### code chunk number 14: BitSeq.Rnw:189-197 ################################################### estimateExpression( "data-c0b0.prob", outFile = "data-c0b0", outputType = "RPKM", trInfoFile = "data.tr", MCMC_burnIn = 200, MCMC_samplesN = 200, MCMC_samplesSave = 100, MCMC_chainsN = 2 ) ################################################### ### code chunk number 15: BitSeq.Rnw:213-222 ################################################### estimateExpressionLegacy( "data-c0b0.prob", outFile = "data-c0b0", outputType = "RPKM", trInfoFile = "data.tr", MCMC_burnIn = 200, MCMC_samplesN = 200, MCMC_samplesSave = 100, MCMC_scaleReduction = 1.1, MCMC_chainsN = 2 ) ################################################### ### code chunk number 16: BitSeq.Rnw:228-229 (eval = FALSE) ################################################### ## loadSamples("data-c0b0.rpkm") ################################################### ### code chunk number 17: BitSeq.Rnw:239-241 ################################################### allConditions = list(c("data-c0b0.rpkm","data-c0b1.rpkm"), c("data-c1b1.rpkm","data-c1b1.rpkm")) ################################################### ### code chunk number 18: BitSeq.Rnw:245-248 ################################################### getMeanVariance(allConditions, outFile = "data.means", log = TRUE ) ################################################### ### code chunk number 19: BitSeq.Rnw:255-259 ################################################### estimateHyperPar( outFile = "data.par", conditions = allConditions, meanFile = "data.means", verbose = TRUE ) ################################################### ### code chunk number 20: BitSeq.Rnw:267-279 ################################################### estimateDE(allConditions, outFile = "data", parFile = "data.par" ) ## ## pretend run with three conditions and normalization constants ## cond3Files = c("data-c2b0.rpkm","data-c2b1.rpkm", "data-c2b2.rpkm") estimateDE(list( allConditions[[1]], allConditions[[2]], cond3Files), outFile="mydata", parFile="mydata.par", norm=c(1.0, 0.999, 1.0017, 0.9998, 1.0, 0.99, 0.97), pretend=TRUE) ################################################### ### code chunk number 21: BitSeq.Rnw:283-287 ################################################### estimateDE(allConditions, outFile = "data", parFile = "data.par", samples = TRUE ) ################################################### ### code chunk number 22: BitSeq.Rnw:300-307 ################################################### res1 <- getExpression("data-c0b0.sam", "ensSelect1.fasta", outPrefix="localDir/data-c0b0", log = TRUE, MCMC_burnIn=200, MCMC_samplesN=200, pretend=TRUE)