BigWigFile-class {rtracklayer} | R Documentation |
BigWig Import and Export
Description
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.
Usage
## S4 method for signature 'BigWigFile,ANY,ANY'
import(con, format, text,
selection = BigWigSelection(which, ...),
which = con,
as = c("GRanges", "RleList", "NumericList"), ...)
import.bw(con, ...)
## S4 method for signature 'ANY,BigWigFile,ANY'
export(object, con, format, ...)
## S4 method for signature 'GenomicRanges,BigWigFile,ANY'
export(object, con, format,
dataFormat = c("auto", "variableStep", "fixedStep",
"bedGraph"), compress = TRUE, fixedSummaries = FALSE)
export.bw(object, con, ...)
Arguments
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. |
as |
Specifies the class of the return object. Default is
|
selection |
A |
which |
A range data structure coercible to |
dataFormat |
Probably best left to “auto”. Exists only for historical reasons. |
compress |
If |
fixedSummaries |
If |
... |
Arguments to pass down to methods to other methods. For
import, the flow eventually reaches the |
Value
A GRanges
(default), RleList
or NumericList
.
GRanges
return ranges with non-zero score values in a score
metadata column. The length of the NumericList
is the same length
as the selection
argument (one list element per range).
The return order in the NumericList
matches the order of the
BigWigSelection
object.
BigWigFile
objects
A BigWigFile
object, an extension of
BiocFile
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 theSeqinfo
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_), as = c("GRangesList", "RleList", "matrix"), ...
: Aggregates the intervals in the file that fall intoranges
, which should be something coercible toGRanges
. The aggregation essentially compresses each sequence to a length ofsize
. The algorithm is specified bytype
; 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 type depends onas
, and can be a GRangesList, RleList or matrix, where the number elements (or rows) is equal to the length ofranges
. Foras="matrix"
, there must be one unique value ofsize
, which is equal to the number of columns in the result. Theas="matrix"
case is the only one that supports asize
greater than the width of the corresponding element inranges
, where values are interpolated to yield the matrix result. 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.If a summary statistic is not available / cannot be computed for a given range a warning is thrown and the defaultValue
NA_real_
is returned.
When accessing remote data, the UCSC library caches data in the
‘/tmp/udcCache’ directory. To clean the cache, call
cleanBigWigCache(maxDays)
, where any files older than
maxDays
days old will be deleted.
BigWigFileList
objects
A BigWigFileList()
provides a convenient way of managing a list
of BigWigFile
instances.
Author(s)
Michael Lawrence
See Also
wigToBigWig
for converting a WIG file to BigWig.
Examples
if (.Platform$OS.type != "windows") {
test_path <- system.file("tests", package = "rtracklayer")
test_bw <- file.path(test_path, "test.bw")
## GRanges
## Returns ranges with non-zero scores.
gr <- import(test_bw)
gr
which <- GRanges(c("chr2", "chr2"), IRanges(c(1, 300), c(400, 1000)))
import(test_bw, which = which)
## RleList
## Scores returned as an RleList is equivalent to the coverage.
## Best option when 'which' or 'selection' contain many small ranges.
mini <- narrow(unlist(tile(which, 50)), 2)
rle <- import(test_bw, which = mini, as = "RleList")
rle
## NumericList
## The 'which' is stored as metadata:
track <- import(test_bw, which = which, as = "NumericList")
metadata(track)
## 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)
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 = seqlengths(bwf) / 10) # 10X reduction
summary(bwf, type = "min") # min instead of mean
summary(bwf, track, size = 10, as = "matrix") # each feature 10 windows
}