### R code from vignette source 'methylationAnalysis.Rnw' ################################################### ### code chunk number 1: <style-Sweave ################################################### BiocStyle::latex() ################################################### ### code chunk number 2: methylationAnalysis.Rnw:28-29 ################################################### library(segmentSeq) ################################################### ### code chunk number 3: methylationAnalysis.Rnw:34-41 ################################################### if(require("parallel")) { numCores <- min(8, detectCores()) cl <- makeCluster(numCores) } else { cl <- NULL } ################################################### ### code chunk number 4: methylationAnalysis.Rnw:46-55 ################################################### datadir <- system.file("extdata", package = "segmentSeq") files <- c("short_18B_C24_C24_trim.fastq_CG_methCalls", "short_Sample_17A_trimmed.fastq_CG_methCalls", "short_13_C24_col_trim.fastq_CG_methCalls", "short_Sample_28_trimmed.fastq_CG_methCalls") mD <- readMeths(files = files, dir = datadir, libnames = c("A1", "A2", "B1", "B2"), replicates = c("A","A","B","B"), nonconversion = c(0.004777, 0.005903, 0.016514, 0.006134)) ################################################### ### code chunk number 5: methDist ################################################### par(mfrow = c(2,1)) dists <- plotMethDistribution(mD, main = "Distributions of methylation", chr = "Chr1") plotMethDistribution(mD, subtract = rowMeans(sapply(dists, function(x) x[,2])), main = "Differences between distributions", chr = "Chr1") ################################################### ### code chunk number 6: figMethDist ################################################### par(mfrow = c(2,1)) dists <- plotMethDistribution(mD, main = "Distributions of methylation", chr = "Chr1") plotMethDistribution(mD, subtract = rowMeans(sapply(dists, function(x) x[,2])), main = "Differences between distributions", chr = "Chr1") ################################################### ### code chunk number 7: methylationAnalysis.Rnw:81-82 ################################################### sD <- processAD(mD, gap = 300, squeeze = 10, filterProp = 0.05, verbose = TRUE, strandSplit = TRUE, cl = cl) ################################################### ### code chunk number 8: methylationAnalysis.Rnw:91-93 ################################################### hS <- heuristicSeg(sD = sD, aD = mD, prop = 0.2, cl = cl, gap = 100, getLikes = FALSE) hS ################################################### ### code chunk number 9: methylationAnalysis.Rnw:98-99 ################################################### hS <- mergeMethSegs(hS, mD, gap = 5000, cl = cl) ################################################### ### code chunk number 10: methylationAnalysis.Rnw:104-105 ################################################### hSL <- lociLikelihoods(hS, mD, cl = cl) ################################################### ### code chunk number 11: plotMeth ################################################### plotMeth(mD, hSL, chr = "Chr1", limits = c(1, 50000), cap = 10) ################################################### ### code chunk number 12: figplotMeth ################################################### plotMeth(mD, hSL, chr = "Chr1", limits = c(1, 50000), cap = 10) ################################################### ### code chunk number 13: methylationAnalysis.Rnw:144-145 ################################################### groups(hSL) <- list(NDE = c(1,1,1,1), DE = c("A", "A", "B", "B")) ################################################### ### code chunk number 14: methylationAnalysis.Rnw:149-150 ################################################### hSL <- methObservables(hSL) ################################################### ### code chunk number 15: methylationAnalysis.Rnw:154-155 ################################################### densityFunction(hSL) <- bbNCDist ################################################### ### code chunk number 16: methylationAnalysis.Rnw:159-160 ################################################### hSL <- getPriors(hSL, cl = cl) ################################################### ### code chunk number 17: methylationAnalysis.Rnw:165-166 ################################################### hSL <- getLikelihoods(hSL, cl = cl) ################################################### ### code chunk number 18: methylationAnalysis.Rnw:171-172 ################################################### topCounts(hSL, "DE") ################################################### ### code chunk number 19: <stopCluster ################################################### if(!is.null(cl)) stopCluster(cl) ################################################### ### code chunk number 20: methylationAnalysis.Rnw:185-186 ################################################### sessionInfo()