The following introduces several utilities useful for ChIP-Seq data. They are not part of the actual workflow.
Rle object stores coverage information
library(rtracklayer)
library(GenomicRanges)
library(Rsamtools)
library(GenomicAlignments)
outpaths <- subsetWF(args, slot = "output", subset = 1, index = 1)
aligns <- readGAlignments(outpaths[1])
cov <- coverage(aligns)
cov
Resizing aligned reads
trim(resize(as(aligns, "GRanges"), width = 200))
Naive peak calling
islands <- slice(cov, lower = 15)
islands[[1]]
Plot coverage for defined region
library(ggbio)
myloc <- c("Chr1", 1, 1e+05)
ga <- readGAlignments(outpaths[1], use.names = TRUE, param = ScanBamParam(which = GRanges(myloc[1],
IRanges(as.numeric(myloc[2]), as.numeric(myloc[3])))))
autoplot(ga, aes(color = strand, fill = strand), facets = strand ~
seqnames, stat = "coverage")

