-
Notifications
You must be signed in to change notification settings - Fork 18
/
Copy pathranges-overlaps.Rd
133 lines (115 loc) · 4.59 KB
/
ranges-overlaps.Rd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/ranges-overlap-find.R,
% R/ranges-overlap-groups.R
\name{find_overlaps}
\alias{find_overlaps}
\alias{find_overlaps.IntegerRanges}
\alias{find_overlaps.GenomicRanges}
\alias{find_overlaps_within}
\alias{find_overlaps_within.IntegerRanges}
\alias{find_overlaps_within.GenomicRanges}
\alias{find_overlaps_directed}
\alias{find_overlaps_directed.GenomicRanges}
\alias{find_overlaps_within_directed}
\alias{find_overlaps_within_directed.GenomicRanges}
\alias{group_by_overlaps}
\alias{group_by_overlaps.IntegerRanges}
\alias{group_by_overlaps.GenomicRanges}
\title{Find overlap between two Ranges}
\usage{
find_overlaps(x, y, maxgap, minoverlap, suffix = c(".x", ".y"))
\method{find_overlaps}{IntegerRanges}(x, y, maxgap = -1L, minoverlap = 0L, suffix = c(".x", ".y"))
\method{find_overlaps}{GenomicRanges}(x, y, maxgap = -1L, minoverlap = 0L, suffix = c(".x", ".y"))
find_overlaps_within(x, y, maxgap, minoverlap, suffix = c(".x", ".y"))
\method{find_overlaps_within}{IntegerRanges}(
x,
y,
maxgap = -1L,
minoverlap = 0L,
suffix = c(".x", ".y")
)
\method{find_overlaps_within}{GenomicRanges}(
x,
y,
maxgap = -1L,
minoverlap = 0L,
suffix = c(".x", ".y")
)
find_overlaps_directed(x, y, maxgap, minoverlap, suffix = c(".x", ".y"))
\method{find_overlaps_directed}{GenomicRanges}(
x,
y,
maxgap = -1L,
minoverlap = 0L,
suffix = c(".x", ".y")
)
find_overlaps_within_directed(x, y, maxgap, minoverlap, suffix = c(".x", ".y"))
\method{find_overlaps_within_directed}{GenomicRanges}(x, y, maxgap, minoverlap, suffix = c(".x", ".y"))
group_by_overlaps(x, y, maxgap, minoverlap)
\method{group_by_overlaps}{IntegerRanges}(x, y, maxgap = -1L, minoverlap = 0L)
\method{group_by_overlaps}{GenomicRanges}(x, y, maxgap = -1L, minoverlap = 0L)
}
\arguments{
\item{x, y}{Objects representing ranges}
\item{maxgap, minoverlap}{The maximimum gap between intervals as an integer
greater than or equal to negative one. The minimum amount of overlap between intervals
as an integer greater than zero, accounting for the maximum gap.}
\item{suffix}{A character vector of length two used to identify metadata columns
coming from x and y.}
}
\value{
A Ranges object with rows corresponding to the
ranges in x that overlap y. In the case of \code{group_by_overlaps()}, returns
a GroupedRanges object, grouped by the number of overlaps
of ranges in x that overlap y (stored in a column called query).
}
\description{
Find overlap between two Ranges
}
\details{
\code{find_overlaps()} will search for any overlaps between ranges
x and y and return a Ranges object of length equal to the number of times x
overlaps y. This Ranges object will have additional metadata columns
corresponding to the metadata columns in y. \code{find_overlaps_within()} is
the same but will only search for overlaps within y. For GRanges objects strand is
ignored, unless \code{find_overlaps_directed()} is used. If the Ranges objects have no
metadata, one could use \code{group_by_overlaps()} to be able to
identify the index of the input Range x that overlaps a Range in y.
Alternatively,
\code{pair_overlaps()} could be used to place the x ranges next to the range
in y they overlap.
}
\examples{
query <- data.frame(start = c(5,10, 15,20), width = 5, gc = runif(4)) \%>\%
as_iranges()
subject <- data.frame(start = 2:6, width = 3:7, label = letters[1:5]) \%>\%
as_iranges()
find_overlaps(query, subject)
find_overlaps(query, subject, minoverlap = 5)
find_overlaps_within(query, subject) # same result as minoverlap
find_overlaps(query, subject, maxgap = 1)
# -- GRanges objects, strand is ignored by default
query <- data.frame(seqnames = "chr1",
start = c(11,101),
end = c(21, 200),
name = c("a1", "a2"),
strand = c("+", "-"),
score = c(1,2)) \%>\%
as_granges()
subject <- data.frame(seqnames = "chr1",
strand = c("+", "-", "+", "-"),
start = c(21,91,101,201),
end = c(30,101,110,210),
name = paste0("b", 1:4),
score = 1:4) \%>\%
as_granges()
# ignores strandedness
find_overlaps(query, subject, suffix = c(".query", ".subject"))
find_overlaps(query, subject, suffix = c(".query", ".subject"), minoverlap = 2)
# adding directed prefix includes strand
find_overlaps_directed(query, subject, suffix = c(".query", ".subject"))
}
\seealso{
\code{IRanges::\link[IRanges:findOverlaps-methods]{findOverlaps()}},
\code{GenomicRanges::\link[GenomicRanges:findOverlaps-methods]{findOverlaps()}}
}