setops-methods {GenomicRanges} | R Documentation |
Performs set operations on GRanges and GRangesList objects.
NOTE: The punion
, pintersect
,
psetdiff
, and pgap
generic
functions and methods for Ranges objects are defined and
documented in the IRanges package.
## Set operations ## S4 method for signature 'GRanges,GRanges' union(x, y, ignore.strand=FALSE, ...) ## S4 method for signature 'GRanges,GRanges' intersect(x, y, ignore.strand=FALSE, ...) ## S4 method for signature 'GRanges,GRanges' setdiff(x, y, ignore.strand=FALSE, ...) ## Parallel set operations ## S4 method for signature 'GRanges,GRanges' punion(x, y, fill.gap=FALSE, ignore.strand=FALSE, ...) ## S4 method for signature 'GRanges,GRanges' pintersect(x, y, drop.nohit.ranges=FALSE, ignore.strand=FALSE, strict.strand=FALSE) ## S4 method for signature 'GRanges,GRanges' psetdiff(x, y, ignore.strand=FALSE, ...)
x, y |
For For For For In addition, for the "parallel" operations, |
fill.gap |
Logical indicating whether or not to force a union by using the rule
|
ignore.strand |
For set operations: If set to TRUE, then the strand of For parallel set operations: If set to TRUE, the strand information is
ignored in the computation and the result has the strand information of
|
drop.nohit.ranges |
If TRUE then elements in If FALSE (the default) then nothing is removed and a |
strict.strand |
If set to FALSE (the default), features on the |
... |
Further arguments to be passed to or from other methods. |
The pintersect
methods involving GRanges and/or
GRangesList objects use the triplet (sequence name, range, strand)
to determine the element by element intersection of features, where a
strand value of "*"
is treated as occurring on both the "+"
and "-"
strand (unless strict.strand
is set to TRUE, in
which case the strand of intersecting elements must be strictly the same).
The psetdiff
methods involving GRanges and/or
GRangesList objects use the triplet (sequence name, range,
strand) to determine the element by element set difference of features,
where a strand value of "*"
is treated as occurring on both the
"+"
and "-"
strand.
For union
, intersect
, setdiff
, and pgap
: a
GRanges object.
For punion
and pintersect
: when x
or y
is
not a GRanges object, an object of the same class as this
non-GRanges object. Otherwise, a GRanges object.
For psetdiff
: either a GRanges object when both x
and y
are GRanges objects, or a GRangesList object
when y
is a GRangesList object.
P. Aboyoun and H. Pagès
setops-methods, GRanges-class, GRangesList-class, findOverlaps-methods
## --------------------------------------------------------------------- ## A. SET OPERATIONS ## --------------------------------------------------------------------- x <- GRanges("chr1", IRanges(c(2, 9) , c(7, 19)), strand=c("+", "-")) y <- GRanges("chr1", IRanges(5, 10), strand="-") union(x, y) union(x, y, ignore.strand=TRUE) intersect(x, y) intersect(x, y, ignore.strand=TRUE) setdiff(x, y) setdiff(x, y, ignore.strand=TRUE) ## --------------------------------------------------------------------- ## B. PARALLEL SET OPERATIONS ## --------------------------------------------------------------------- punion(x, shift(x, 6)) ## Not run: punion(x, shift(x, 7)) # will fail ## End(Not run) punion(x, shift(x, 7), fill.gap=TRUE) pintersect(x, shift(x, 6)) pintersect(x, shift(x, 7)) psetdiff(x, shift(x, 7)) ## --------------------------------------------------------------------- ## C. MORE EXAMPLES ## --------------------------------------------------------------------- ## GRanges object: gr <- GRanges(seqnames=c("chr2", "chr1", "chr1"), ranges=IRanges(1:3, width = 12), strand=Rle(strand(c("-", "*", "-")))) ## GRangesList object gr1 <- GRanges(seqnames="chr2", ranges=IRanges(3, 6)) 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("-", "-")) grlist <- GRangesList(gr1=gr1, gr2=gr2, gr3=gr3) ## Parallel intersection of a GRanges and a GRangesList object pintersect(gr, grlist) pintersect(grlist, gr) ## For a fast 'mendoapply(intersect, grlist, as(gr, "GRangesList"))' ## call pintersect() with 'strict.strand=TRUE' and call reduce() on ## the result with 'drop.empty.ranges=TRUE': reduce(pintersect(grlist, gr, strict.strand=TRUE), drop.empty.ranges=TRUE) ## Parallel set difference of a GRanges and a GRangesList object psetdiff(gr, grlist) ## Parallel set difference of two GRangesList objects psetdiff(grlist, shift(grlist, 3))