BigWigFile-class {rtracklayer} | R Documentation |
These functions support the import and export of the UCSC BigWig format, a compressed, binary form of WIG/BEDGraph with a spatial index and precomputed summaries. These functions do not work on Windows.
## S4 method for signature 'BigWigFile,ANY,ANY' import(con, format, text, selection = BigWigSelection(which, ...), which = con, asRangedData = FALSE, ...) import.bw(con, ...) ## S4 method for signature 'ANY,BigWigFile,ANY' export(object, con, format, ...) ## S4 method for signature 'RangedData,BigWigFile,ANY' export(object, con, format, dataFormat = c("auto", "variableStep", "fixedStep", "bedGraph"), compress = TRUE) export.bw(object, con, ...)
con |
A path, URL or |
object |
The object to export, should be an |
format |
If not missing, should be “bigWig” or “bw” (case insensitive). |
text |
Not supported. |
asRangedData |
If |
selection |
A |
which |
A range data structure coercible to |
dataFormat |
Probably best left to “auto”. Exists only for historical reasons. |
compress |
If |
... |
Arguments to pass down to methods to other methods. For
import, the flow eventually reaches the |
A GRanges
(or RangedData
if asRangedData
is
TRUE
), with the score values in the score
metadata column,
which is accessible via the score
function.
BigWigFile
objectsA BigWigFile
object, an extension of
RTLFile
is a reference to a BigWig file. To cast
a path, URL or connection to a BigWigFile
, pass it to the
BigWigFile
constructor.
BigWig files are more complex than most track files, and there are a
number of methods on BigWigFile
for accessing the additional
information:
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.
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
wigToBigWig
for converting a WIG file to BigWig.
if (.Platform$OS.type != "windows") { test_path <- system.file("tests", package = "rtracklayer") test_bw <- file.path(test_path, "test.bw") test <- import(test_bw, asRangedData = FALSE) test which <- GRanges(c("chr2", "chr2"), IRanges(c(1, 300), c(400, 1000))) import(test_bw, which = which, asRangedData = FALSE) ## Not run: test_bw_out <- file.path(tempdir(), "test_out.bw") export(test, test_bw_out) ## End(Not run) bwf <- BigWigFile(test_bw) track <- import(bwf, asRangedData = FALSE) seqinfo(bwf) 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 }