| coverage-methods {GenomicRanges} | R Documentation |
coverage methods for GRanges,
GRangesList and GappedAlignments objects.
## S4 method for signature 'GenomicRanges' coverage(x, shift=0L, width=NULL, weight=1L, ...)
x |
A GRanges, GRangesList or GappedAlignments object. |
shift, width, weight, ... |
See |
Here is how optional arguments shift, width and
weight are handled when x is a GRanges object:
shift, weight: can be either a numeric vector
(integers) or a list. If a list, then it should be named by the
sequence levels in x (i.e. by the names of the underlying
sequences), and its elements are passed into the coverage
method for IRanges objects. If a numeric vector,
then it is first recycled to the length of x, then turned
into a list with split(shift, as.factor(seqnames(x))),
and finally the elements of this list are passed into the
coverage method for IRanges objects.
width: can be either NULL or a numeric vector.
If a numeric vector, then it should be named by the sequence
levels in x. If NULL (the default), then it is
replaced with seqlengths(x). Like for shift and
weight, its elements are passed into the coverage
method for IRanges objects (if the element is
NA then NULL is passed instead).
When x is a GRangesList object, coverage(x, ...)
is equivalent to coverage(unlist(x), ...).
When x is a GappedAlignments object, coverage(x, ...)
is equivalent to coverage(as(x, "GRangesList"), ...).
Returns a named RleList object with one element
('integer' Rle) per underlying sequence in x representing how
many times each position in the sequence is covered by the intervals in
x.
P. Aboyoun and H. Pages
coverage,
RleList-class,
GRanges-class,
GRangesList-class,
GappedAlignments-class
## Coverage of a GRanges object:
gr <- GRanges(
seqnames=Rle(c("chr1", "chr2", "chr1", "chr3"), c(1, 3, 2, 4)),
ranges=IRanges(1:10, end=10),
strand=Rle(strand(c("-", "+", "*", "+", "-")), c(1, 2, 2, 3, 2)),
seqlengths=c(chr1=11, chr2=12, chr3=13))
cvg <- coverage(gr)
pcvg <- coverage(gr[strand(gr) == "+"])
mcvg <- coverage(gr[strand(gr) == "-"])
scvg <- coverage(gr[strand(gr) == "*"])
stopifnot(identical(pcvg + mcvg + scvg, cvg))
## Coverage of a GRangesList object:
gr1 <- GRanges(seqnames="chr2",
ranges=IRanges(3, 6),
strand = "+")
gr2 <- GRanges(seqnames=c("chr1", "chr1"),
ranges=IRanges(c(7,13), width=3),
strand=c("+", "-"))
gr3 <- GRanges(seqnames=c("chr1", "chr2"),
ranges=IRanges(c(1, 4), c(3, 9)),
strand=c("-", "-"))
grl <- GRangesList(gr1=gr1, gr2=gr2, gr3=gr3)
stopifnot(identical(coverage(grl), coverage(unlist(grl))))
## Coverage of a GappedAlignments object:
library(Rsamtools) # because file ex1.bam is in this package
galn_file <- system.file("extdata", "ex1.bam", package="Rsamtools")
galn <- readGappedAlignments(galn_file)
stopifnot(identical(coverage(galn), coverage(as(galn, "GRangesList"))))