if (all(requireNamespace("RNAseqData.HNRNPC.bam.chr14", quietly=TRUE) &&
require(GenomicAlignments))) {
## -----------------------------------------------------------------------
## Compute coverage across BAM files.
## -----------------------------------------------------------------------
fls <- ## 8 bam files
RNAseqData.HNRNPC.bam.chr14::RNAseqData.HNRNPC.bam.chr14_BAMFILES
## Regions of interest.
gr <- GRanges("chr14", IRanges(c(62262735, 63121531, 63980327),
width=214700))
## The MAP computes the coverage ...
MAP <- function(range, file, ...) {
requireNamespace("GenomicFiles", quietly=TRUE)
## for coverage(), Rsamtools::ScanBamParam()
param = Rsamtools::ScanBamParam(which=range)
GenomicFiles::coverage(file, param=param)[range]
}
## and REDUCE adds the last and current results.
REDUCE <- function(mapped, ...)
Reduce("+", mapped)
## -----------------------------------------------------------------------
## reduceByRange:
## With no REDUCE, coverage is computed for each range / file combination.
cov1 <- reduceByRange(gr, fls, MAP)
cov1[[1]]
## Each call to coverage() produces an RleList which accumulate on the
## workers. We can use a reducer to combine these lists either iteratively
## or non-iteratively. When iterate = TRUE the current result
## is collapsed with the last resulting in a maximum of 2 RleLists on
## a worker at any given time.
cov2 <- reduceByRange(gr, fls, MAP, REDUCE, iterate=TRUE)
cov2[[1]]
## If memory use is not a concern (or if MAP output is not large) the
## REDUCE function can be applied non-iteratively.
cov3 <- reduceByRange(gr, fls, MAP, REDUCE, iterate=FALSE)
## Results match those obtained with the iterative REDUCE.
cov3[[1]]
## When 'ranges' is a GRangesList, the list elements are sent to the
## workers instead of a single range as in the case of a GRanges.
grl <- GRangesList(gr[1], gr[2:3])
grl
cov4 <- reduceByRange(grl, fls, MAP)
length(cov4) ## length of GRangesList
elementNROWS(cov4) ## number of files
## -----------------------------------------------------------------------
## reduceRanges:
## This function passes the character vector of all file names to MAP.
## MAP must handle each file separately or invoke a method that operates
## on a list of files.
## TODO: example
}
Run the code above in your browser using DataLab