1
\name{coverage-methods}
2
\alias{coverage-methods}
4
\alias{coverage,GenomicRanges-method}
5
\alias{coverage,GRangesList-method}
6
\alias{coverage,GappedAlignments-method}
7
\alias{coverage,GappedAlignmentPairs-method}
10
\title{Coverage of a GRanges, GRangesList, GappedAlignments, or
11
GappedAlignmentPairs object}
14
\code{\link[IRanges]{coverage}} methods for \link{GRanges},
15
\link{GRangesList}, \link{GappedAlignments}, and \link{GappedAlignmentPairs}
20
\S4method{coverage}{GenomicRanges}(x, shift=0L, width=NULL, weight=1L, ...)
21
\S4method{coverage}{GappedAlignments}(x, shift=0L, width=NULL,
22
weight=1L, drop.D.ranges=FALSE, ...)
23
\S4method{coverage}{GappedAlignmentPairs}(x, shift=0L, width=NULL,
24
weight=1L, drop.D.ranges=FALSE, ...)
29
A \link{GRanges}, \link{GRangesList}, \link{GappedAlignments},
30
or \link{GappedAlignmentPairs} object.
32
\item{shift, width, weight, ...}{
33
See \code{\link[IRanges]{coverage}} in the IRanges package for
34
a description of these optional arguments.
37
Whether the coverage calculation should ignore ranges corresponding
38
to D (deletion) in the CIGAR string.
43
Here is how optional arguments \code{shift}, \code{width} and
44
\code{weight} are handled when \code{x} is a \link{GRanges} object:
46
\item \code{shift}, \code{weight}: can be either a numeric vector
47
(integers) or a list. If a list, then it should be named by the
48
sequence levels in \code{x} (i.e. by the names of the underlying
49
sequences), and its elements are passed into the \code{coverage}
50
method for \link[IRanges]{IRanges} objects. If a numeric vector,
51
then it is first recycled to the length of \code{x}, then turned
52
into a list with \code{split(shift, as.factor(seqnames(x)))},
53
and finally the elements of this list are passed into the
54
\code{coverage} method for \link[IRanges]{IRanges} objects.
55
Finally, if \code{x} is a \link{GRanges} object, then \code{weight}
56
can also be a single string naming a metadata column to be used
59
\item \code{width}: can be either \code{NULL} or a numeric vector.
60
If a numeric vector, then it should be named by the sequence
61
levels in \code{x}. If \code{NULL} (the default), then it is
62
replaced with \code{seqlengths(x)}. Like for \code{shift} and
63
\code{weight}, its elements are passed into the \code{coverage}
64
method for \link[IRanges]{IRanges} objects (if the element is
65
\code{NA} then \code{NULL} is passed instead).
68
When \code{x} is a \link{GRangesList} object, \code{coverage(x, ...)}
69
is equivalent to \code{coverage(unlist(x), ...)}.
71
When \code{x} is a \link{GappedAlignments} or \link{GappedAlignmentPairs}
72
object, \code{coverage(x, ...)} is equivalent to
73
\code{coverage(as(x, "GRangesList"), ...)}.
77
Returns a named \link[IRanges]{RleList} object with one element
78
('integer' Rle) per underlying sequence in \code{x} representing how
79
many times each position in the sequence is covered by the intervals in
83
\author{P. Aboyoun and H. Pages}
87
\item \code{\link[IRanges]{coverage}}.
88
\item \link[IRanges]{RleList-class}.
89
\item \link{GRanges-class}.
90
\item \link{GRangesList-class}.
91
\item \link{GappedAlignments-class}.
92
\item \link{GappedAlignmentPairs-class}.
97
## Coverage of a GRanges object:
99
seqnames=Rle(c("chr1", "chr2", "chr1", "chr3"), c(1, 3, 2, 4)),
100
ranges=IRanges(1:10, end=10),
101
strand=Rle(strand(c("-", "+", "*", "+", "-")), c(1, 2, 2, 3, 2)),
102
seqlengths=c(chr1=11, chr2=12, chr3=13))
104
pcvg <- coverage(gr[strand(gr) == "+"])
105
mcvg <- coverage(gr[strand(gr) == "-"])
106
scvg <- coverage(gr[strand(gr) == "*"])
107
stopifnot(identical(pcvg + mcvg + scvg, cvg))
109
## Coverage of a GRangesList object:
110
gr1 <- GRanges(seqnames="chr2",
111
ranges=IRanges(3, 6),
113
gr2 <- GRanges(seqnames=c("chr1", "chr1"),
114
ranges=IRanges(c(7,13), width=3),
116
gr3 <- GRanges(seqnames=c("chr1", "chr2"),
117
ranges=IRanges(c(1, 4), c(3, 9)),
119
grl <- GRangesList(gr1=gr1, gr2=gr2, gr3=gr3)
120
stopifnot(identical(coverage(grl), coverage(unlist(grl))))
122
## Coverage of a GappedAlignments or GappedAlignmentPairs object:
123
library(Rsamtools) # because file ex1.bam is in this package
124
ex1_file <- system.file("extdata", "ex1.bam", package="Rsamtools")
125
galn <- readGappedAlignments(ex1_file)
126
stopifnot(identical(coverage(galn), coverage(as(galn, "GRangesList"))))
127
galp <- readGappedAlignmentPairs(ex1_file)
128
stopifnot(identical(coverage(galp), coverage(as(galp, "GRangesList"))))