Commit b920ead9 authored by Safia Saci's avatar Safia Saci
Browse files

Update bayMethList_GD.R

parent 2ddfefcf
#!/usr/bin/env Rscript
suppressWarnings(suppressPackageStartupMessages(library(optparse)))
suppressWarnings(suppressPackageStartupMessages(library(Repitools)))
suppressWarnings(suppressPackageStartupMessages(library(BSgenome.Sscrofa.UCSC.susScr11)))
suppressWarnings(suppressPackageStartupMessages(library(BSgenome.Sscrofa.ENSEMBL.susScr11)))
suppressWarnings(suppressPackageStartupMessages(library(stringr)))
option_list <- list(
......@@ -32,28 +32,32 @@ results <- BayMethList(
results <- determineOffset(results, controlPlot=list(show=FALSE))
results <- empBayes(results,
verbose = TRUE,
ngroups = 100,
ncomp = 1,
maxBins = length(results)%/%10,
method="beta",
ncpu = 1)
verbose = TRUE,
ngroups = 100,
ncomp = 1,
maxBins = length(results)%/%10,
method="beta",
ncpu = 1)
results <- methylEst(results,
verbose = TRUE,
controlCI = list(compute = FALSE,
method = "Wald",
level = 0.95,
nmarg = 300,
ncpu = 1))
verbose = TRUE,
controlCI = list(compute = FALSE,
method = "Wald",
level = 0.95,
nmarg = 300,
ncpu = 1))
methTrack <- windows_genome
ssLengths <- seqlengths(Sscrofa)[paste0("chr", c(1:18, "X", "Y"))]
names(ssLengths) <- str_replace(names(ssLengths), "chr", "")
seqlengths(methTrack) <- ssLengths
methTrack$score <- methEst(results)$mean[, 1]
methTrack$score[is.na(methTrack$score)] <- 0
seqlengths(methTrack) <- seqlengths(Sscrofa)
export(methTrack, opt$output_name, format = "bigWig")
methTrack2 <- windows_genome
methTrack2$score <- results@sampleInterest[, 1]
methTrack2$score[is.na(methTrack2$score)] <- 0
seqlengths(methTrack2) <- seqlengths(Sscrofa)
export(methTrack2, paste0("raw_count_", opt$output_name), format = "bigWig")
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment