### R code from vignette source 'MethylSeekR.Rnw' ################################################### ### code chunk number 1: install (eval = FALSE) ################################################### ## if (!requireNamespace("BiocManager", quietly=TRUE)) ## install.packages("BiocManager") ## BiocManager::install("BSgenome") ################################################### ### code chunk number 2: genomes (eval = FALSE) ################################################### ## library(BSgenome) ## available.genomes() ################################################### ### code chunk number 3: install (eval = FALSE) ################################################### ## if (!requireNamespace("BiocManager", quietly=TRUE)) ## install.packages("BiocManager") ## BiocManager::install("BSgenome.Hsapiens.UCSC.hg18") ################################################### ### code chunk number 4: MethylSeekR.Rnw:101-102 ################################################### library(MethylSeekR) ################################################### ### code chunk number 5: MethylSeekR.Rnw:111-112 ################################################### set.seed(123) ################################################### ### code chunk number 6: MethylSeekR.Rnw:127-129 (eval = FALSE) ################################################### ## system.file("extdata", "Lister2009_imr90_hg18_chr22.tab", ## package="MethylSeekR") ################################################### ### code chunk number 7: MethylSeekR.Rnw:139-142 ################################################### library("BSgenome.Hsapiens.UCSC.hg18") sLengths=seqlengths(Hsapiens) head(sLengths) ################################################### ### code chunk number 8: MethylSeekR.Rnw:158-162 ################################################### methFname <- system.file("extdata", "Lister2009_imr90_hg18_chr22.tab", package="MethylSeekR") meth.gr <- readMethylome(FileName=methFname, seqLengths=sLengths) head(meth.gr) ################################################### ### code chunk number 9: MethylSeekR.Rnw:188-190 (eval = FALSE) ################################################### ## system.file("extdata", "SNVs_hg18_chr22.tab", ## package="MethylSeekR") ################################################### ### code chunk number 10: MethylSeekR.Rnw:197-201 ################################################### snpFname <- system.file("extdata", "SNVs_hg18_chr22.tab", package="MethylSeekR") snps.gr <- readSNPTable(FileName=snpFname, seqLengths=sLengths) head(snps.gr) ################################################### ### code chunk number 11: MethylSeekR.Rnw:211-212 ################################################### meth.gr <- removeSNPs(meth.gr, snps.gr) ################################################### ### code chunk number 12: MethylSeekR.Rnw:235-237 ################################################### plotAlphaDistributionOneChr(m=meth.gr, chr.sel="chr22", num.cores=1) ################################################### ### code chunk number 13: MethylSeekR.Rnw:255-257 (eval = FALSE) ################################################### ## library(parallel) ## detectCores() ################################################### ### code chunk number 14: MethylSeekR.Rnw:269-272 ################################################### PMDsegments.gr <- segmentPMDs(m=meth.gr, chr.sel="chr22", seqLengths=sLengths, num.cores=1) head(PMDsegments.gr) ################################################### ### code chunk number 15: MethylSeekR.Rnw:291-294 ################################################### plotAlphaDistributionOneChr(m=subsetByOverlaps(meth.gr, PMDsegments.gr[values(PMDsegments.gr)$type=="notPMD"]), chr.sel="chr22", num.cores=1) ################################################### ### code chunk number 16: myfig1 ################################################### plotPMDSegmentation(m=meth.gr, segs=PMDsegments.gr) ################################################### ### code chunk number 17: MethylSeekR.Rnw:332-334 (eval = FALSE) ################################################### ## savePMDSegments(PMDs=PMDsegments.gr, ## GRangesFilename="PMDs.gr.rds", TableFilename="PMDs.tab") ################################################### ### code chunk number 18: MethylSeekR.Rnw:375-381 ################################################### library(rtracklayer) session <- browserSession() genome(session) <- "hg18" query <- ucscTableQuery(session, "cpgIslandExt") CpGislands.gr <- track(query) genome(CpGislands.gr) <- NA ################################################### ### code chunk number 19: MethylSeekR.Rnw:391-393 ################################################### CpGislands.gr <- suppressWarnings(resize(CpGislands.gr, 5000, fix="center")) ################################################### ### code chunk number 20: MethylSeekR.Rnw:398-401 ################################################### stats <- calculateFDRs(m=meth.gr, CGIs=CpGislands.gr, PMDs=PMDsegments.gr, num.cores=1) stats ################################################### ### code chunk number 21: MethylSeekR.Rnw:424-429 ################################################### FDR.cutoff <- 5 m.sel <- 0.5 n.sel=as.integer(names(stats$FDRs[as.character(m.sel), ] [stats$FDRs[as.character(m.sel), ]