coverage-methods {GenomicRanges}R Documentation

Coverage of a GRanges, GRangesList, GappedAlignments, or GappedAlignmentPairs object

Description

coverage methods for GRanges, GRangesList, GappedAlignments, and GappedAlignmentPairs objects.

Usage

## S4 method for signature 'GenomicRanges'
coverage(x, shift=0L, width=NULL, weight=1L, ...)
## S4 method for signature 'GappedAlignments'
coverage(x, shift=0L, width=NULL,
         weight=1L, drop.D.ranges=FALSE, ...)
## S4 method for signature 'GappedAlignmentPairs'
coverage(x, shift=0L, width=NULL,
         weight=1L, drop.D.ranges=FALSE, ...)

Arguments

x

A GRanges, GRangesList, GappedAlignments, or GappedAlignmentPairs object.

shift, width, weight, ...

See coverage in the IRanges package for a description of these optional arguments.

drop.D.ranges

Whether the coverage calculation should ignore ranges corresponding to D (deletion) in the CIGAR string.

Details

Here is how optional arguments shift, width and weight are handled when x is a GRanges object:

When x is a GRangesList object, coverage(x, ...) is equivalent to coverage(unlist(x), ...).

When x is a GappedAlignments or GappedAlignmentPairs object, coverage(x, ...) is equivalent to coverage(as(x, "GRangesList"), ...).

Value

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.

Author(s)

P. Aboyoun and H. Pages

See Also

Examples

## 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 or GappedAlignmentPairs object:
library(Rsamtools)  # because file ex1.bam is in this package
ex1_file <- system.file("extdata", "ex1.bam", package="Rsamtools")
galn <- readGappedAlignments(ex1_file)
stopifnot(identical(coverage(galn), coverage(as(galn, "GRangesList"))))
galp <- readGappedAlignmentPairs(ex1_file)
stopifnot(identical(coverage(galp), coverage(as(galp, "GRangesList"))))

[Package GenomicRanges version 1.12.1 Index]