1
\name{GenomicRanges-comparison}
3
\alias{GenomicRanges-comparison}
5
\alias{compare,GenomicRanges,GenomicRanges-method}
7
\alias{==,GenomicRanges,GenomicRanges-method}
8
\alias{<=,GenomicRanges,GenomicRanges-method}
10
\alias{duplicated,GenomicRanges-method}
11
\alias{duplicated.GenomicRanges}
13
\alias{match,GenomicRanges,GenomicRanges-method}
14
\alias{\%in\%,GenomicRanges,GenomicRanges-method}
16
\alias{findMatches,GenomicRanges,GenomicRanges-method}
17
\alias{countMatches,GenomicRanges,GenomicRanges-method}
19
\alias{order,GenomicRanges-method}
20
\alias{rank,GenomicRanges-method}
23
\title{Comparing and ordering genomic ranges}
26
Methods for comparing and ordering the elements in one or more
27
\link{GenomicRanges} objects.
31
## Element-wise (aka "parallel") comparison of 2 GenomicRanges objects
32
## -------------------------------------------------------------------
34
\S4method{==}{GenomicRanges,GenomicRanges}(e1, e2)
36
\S4method{<=}{GenomicRanges,GenomicRanges}(e1, e2)
41
\S4method{duplicated}{GenomicRanges}(x, incomparables=FALSE, fromLast=FALSE,
42
method=c("auto", "quick", "hash"))
47
\S4method{match}{GenomicRanges,GenomicRanges}(x, table, nomatch=NA_integer_, incomparables=NULL,
48
method=c("auto", "quick", "hash"), ignore.strand=FALSE, match.if.overlap=FALSE)
50
## order() and related methods
51
## ----------------------------
53
\S4method{order}{GenomicRanges}(..., na.last=TRUE, decreasing=FALSE)
55
\S4method{rank}{GenomicRanges}(x, na.last=TRUE,
56
ties.method=c("average", "first", "random", "max", "min"))
58
## Generalized element-wise (aka "parallel") comparison of 2 GenomicRanges
60
## ------------------------------------------------------------------------
62
\S4method{compare}{GenomicRanges,GenomicRanges}(x, y)
66
\item{e1, e2, x, table, y}{
67
\link{GenomicRanges} objects.
72
\item{fromLast, method, nomatch}{
73
See \code{?`\link[IRanges]{Ranges-comparison}`} in the IRanges
74
package for a description of these arguments.
77
Whether or not the strand should be ignored when comparing 2 genomic
80
\item{match.if.overlap}{
81
You temporarily need to explicitly provide this argument otherwise
82
\code{match()} will issue a warning.
84
Starting with BioC 2.12, the default behavior of \code{match()} on
85
\link{GenomicRanges} objects has changed to use \emph{equality} instead
86
of \emph{overlap} for comparing elements between \link{GenomicRanges}
87
objects \code{x} and \code{table}. Now \code{x[i]} and \code{table[j]} are
88
considered to match when they are equal (i.e. \code{x[i] == table[j]}),
89
instead of when they overlap. This new behavior is consistent with
90
\code{base::\link[base]{match}()}.
92
If you need the new behavior, use \code{match.if.overlap=FALSE}.
94
If you need the old behavior, you can either do:
95
\code{findOverlaps(x, table, select="first")}
96
which is the recommended way. Alternatively, you can call \code{match()}
97
with \code{match.if.overlap=TRUE}.
100
Additional \link{GenomicRanges} objects used for breaking ties.
106
\code{TRUE} or \code{FALSE}.
109
A character string specifying how ties are treated. Only \code{"first"}
110
is supported for now.
115
Two elements of a \link{GenomicRanges} object (i.e. two genomic ranges) are
116
considered equal iff they are on the same underlying sequence and strand,
117
and have the same start and width. \code{duplicated()} and \code{unique()}
118
on a \link{GenomicRanges} object are conforming to this.
120
The "natural order" for the elements of a \link{GenomicRanges} object is to
121
order them (a) first by sequence level, (b) then by strand, (c) then by
122
start, (d) and finally by width.
123
This way, the space of genomic ranges is totally ordered.
124
Note that the \code{reduce} method for \link{GenomicRanges} uses this
125
"natural order" implicitly. Also, note that, because we already do (c)
126
and (d) for regular ranges (see \code{?`\link[IRanges]{Ranges-comparison}`}),
127
genomic ranges that belong to the same underlying sequence and strand are
128
ordered like regular ranges.
129
\code{order()}, \code{sort()}, and \code{rank()} on a \link{GenomicRanges}
130
object are using this "natural order".
132
Also the \code{==}, \code{!=}, \code{<=}, \code{>=}, \code{<} and \code{>}
133
operators between 2 \link{GenomicRanges} objects are using this "natural
141
\item The \link{GenomicRanges} class.
142
\item \link[IRanges]{Ranges-comparison} in the IRanges
143
package for comparing and ordering genomic ranges.
144
\item \link[GenomicRanges]{intra-range-methods} and
145
\link[GenomicRanges]{inter-range-methods} for
146
intra and inter range transformations.
147
\item \link[GenomicRanges]{setops-methods} for set operations on
148
\link{GenomicRanges} objects.
149
\item \link[GenomicRanges]{findOverlaps-methods} for finding
150
overlapping genomic ranges.
156
seqnames=Rle(c("chr1", "chr2", "chr1", "chr3"), c(1, 3, 2, 4)),
157
ranges=IRanges(c(1:9,7L), end=10),
158
strand=Rle(strand(c("-", "+", "*", "+", "-")), c(1, 2, 2, 3, 2)),
159
seqlengths=c(chr1=11, chr2=12, chr3=13))
160
gr <- c(gr0, gr0[7:3])
162
## ---------------------------------------------------------------------
163
## A. ELEMENT-WISE (AKA "PARALLEL") COMPARISON OF 2 GenomicRanges OBJECTS
164
## ---------------------------------------------------------------------
165
gr[2] == gr[2] # TRUE
166
gr[2] == gr[5] # FALSE
170
## ---------------------------------------------------------------------
171
## B. duplicated(), unique()
172
## ---------------------------------------------------------------------
176
## ---------------------------------------------------------------------
178
## ---------------------------------------------------------------------
180
match(gr, table) # Warning! The warning will be removed in BioC 2.13.
182
## In the meantime, specify 'match.if.overlap=FALSE' to suppress the warning:
183
match(gr, table, match.if.overlap=FALSE)
184
match(gr, table, ignore.strand=TRUE, match.if.overlap=FALSE)
186
gr \%in\% table # Warning! The warning will be removed in BioC 2.13.
187
## In the meantime, use suppressWarnings() to suppress the warning:
188
suppressWarnings(gr \%in\% table)
190
## ---------------------------------------------------------------------
191
## D. findMatches(), countMatches()
192
## ---------------------------------------------------------------------
193
findMatches(gr, table)
194
countMatches(gr, table)
196
findMatches(gr, table, ignore.strand=TRUE)
197
countMatches(gr, table, ignore.strand=TRUE)
199
gr_levels <- unique(gr)
200
countMatches(gr_levels, gr)
202
## ---------------------------------------------------------------------
203
## E. order() AND RELATED METHODS
204
## ---------------------------------------------------------------------
209
## ---------------------------------------------------------------------
210
## F. GENERALIZED ELEMENT-WISE COMPARISON OF 2 GenomicRanges OBJECTS
211
## ---------------------------------------------------------------------
212
gr2 <- GRanges(c(rep("chr1", 12), "chr2"), IRanges(c(1:11, 6:7), width=3))
213
strand(gr2)[12] <- "+"
214
gr3 <- GRanges("chr1", IRanges(5, 9))
217
rangeComparisonCodeToLetter(compare(gr2, gr3))