BigWigFile-class {rtracklayer} | R Documentation |
A BigWigFile
object is a reference to a BigWig file. It exists
to support methods with behavior particular to BigWig files.
In the code snippets below, x
represents a BigWigFile
object.
seqinfo(x)
:
Gets the Seqinfo
object
indicating the lengths of the sequences for the intervals in the
file. No circularity or genome information is available.
import.bw(con, selection = BigWigSelection(ranges, ...),
ranges = con, ...)
: Imports the intervals from a big wig file
con
, according to selection
, a
RangedSelection object indicating the intervals to
retrieve from a bigWig file. Supported types of con
include a BigWigFile
and a file name. Note that this
retrieval is very efficient, due to the indexing of the bigWig
format.
summary(ranges = as(seqinfo(object), "GenomicRanges"), size
= 1L, type = c("mean", "min", "max", "coverage", "sd"),
defaultValue = NA_real_)
: Aggregates the intervals in the file
that fall into ranges
, which should be something
coercible to GRanges
. The aggregation essentially
compresses each sequence to a length of size
. The
algorithm is specified by type
; available algorithms
include the mean, min, max, coverage (percent sequence covered
by at least one feature), and standard deviation. When a window
contains no features, defaultValue
is assumed. The result
is an RleList
, with an
element for each element in ranges
. The
driving use case for this is visualization of coverage when the
screen space is small compared to the viewed portion of the
sequence. The operation is very fast, as it leverages cached
multi-level summaries present in every BigWig file.
Michael Lawrence
import.bw
and export.bw
for reading and
writing BigWig files, respectively.
bwf <- BigWigFile(system.file("tests", "test.bw", package = "rtracklayer")) seqinfo(bwf) track <- import.bw(bwf, asRangedData = FALSE) summary(bwf) # for each sequence, average all values into one summary(bwf, range(head(track))) # just average the first few features summary(bwf, size = GenomicRanges::seqlengths(bwf) / 10) # 10X reduction summary(bwf, type = "min") # min instead of mean