consolidateWindows {csaw} | R Documentation |
Consolidate windows of different sizes into a common set of clusters.
consolidateWindows(data.list, equiweight=TRUE, merge.args=list(), region=NULL, overlap.args=list(), sign.list=NULL) # Deprecated, use consolidateTests or consolidateOverlaps instead! consolidateSizes(data.list, result.list, equiweight=TRUE, merge.args=list(), region=NULL, overlap.args=list(), ...)
data.list |
A list of RangedSummarizedExperiment and/or GRanges objects. |
result.list |
A list of data.frames containing the DB test results for each entry of |
equiweight |
A logical scalar indicating whether equal weighting from each window size should be enforced. |
merge.args |
A named list of parameters to pass to |
region |
A GRanges object specifying regions of interest for overlapping with windows. |
overlap.args |
A named list of parameters to pass to |
sign.list |
A list of logical vectors, where each vector corresponds to and is of the same length as an object in |
... |
Further arguments to pass to |
This function merges windows of differing sizes, generated by multiple calls to windowCounts
using different width
values.
If region=NULL
, windows of all sizes are clustered together through mergeWindows
, using a default tolerance of 100 bp.
Otherwise, clusters are defined by overlapping windows of all sizes to region
using findOverlaps
, where each entry of region
defines a cluster.
The aim is to consolidate statistics across window sizes using consolidateTests
or consolidateOverlaps
.
This will provide comprehensive detection of DB at a range of spatial resolutions.
However, this requires some care to balance the contribution of analyses with different window sizes to the combined p-value.
This is achieved by weighting each window and using the weights in downstream functions.
By default, we set equiweight=TRUE
to offset the differences in the number of windows of each size per cluster.
The weight of each window is inversely proportional to the number of windows of the same size within a given cluster.
If region!=NULL
, the weight of each overlap is inversely proportional to the number of overlaps involving windows of the same size to the given region.
This ensures that the combined p-value is not fully determined by numerous small windows within a cluster.
The signs of the log-fold change can be specified using sign.list
.
This should be a list of logical vectors indicating whether the log-fold changes for each set of windows of data.list
are positive.
We do not recommend using this for routine analyses - see ?mergeWindows
for details.
If region=NULL
, a named list is returned containing:
id
:A list of integer vectors, where each vector corresponds to an object in data.list
.
The entries of the vector specify the cluster to which each window of that object is assigned.
weight
:A list of numeric vectors, where each vector corresponds to a vector in id
.
The entries of the vector specify the weight of the corresponding windows.
region
:A GRanges object containing the genomic coordinates of the clusters of merged windows.
If region
is not NULL
, a named list is returned containing:
olap
:A list of Hits object, where each object corresponds to an object in data.list
.
Each Hits object contains the overlaps between the input regions
(query) and the windows in the corresponding object in data.list
(subject).
weight
:A list of numeric vectors, where each vector corresponds to a Hits object in olap
.
The entries of the vector specify the per-overlap weight of the corresponding windows.
Aaron Lun
windowCounts
,
mergeWindows
,
findOverlaps
,
consolidateTests
,
consolidateOverlaps
bamFiles <- system.file("exdata", c("rep1.bam", "rep2.bam"), package="csaw") low <- windowCounts(bamFiles, width=1, filter=1) med <- windowCounts(bamFiles, width=100, filter=1) high <- windowCounts(bamFiles, width=500, filter=1) # Consolidating. cons <- consolidateWindows(list(low, med, high), merge.args=list(tol=100)) str(cons$id) str(cons$weights) cons$region # Trying with a custom region. of.interest <- GRanges(c('chrA', 'chrA', 'chrB', 'chrC'), IRanges(c(1, 500, 100, 1000), c(200, 1000, 700, 1500))) cons <- consolidateWindows(list(low, med, high), region=of.interest) str(cons$id) str(cons$weights)