Description Usage Arguments Details Author(s) See Also Examples
This man page documents intra range transformations of a GenomicRanges object (i.e. of an object that belongs to the GenomicRanges class or one of its subclasses, this includes for example GRanges objects), or a GRangesList object.
See ?`intra-range-methods`
and
?`inter-range-methods`
in the IRanges
package for a quick introduction to intra range and inter
range transformations.
Intra range methods for GAlignments and GAlignmentsList objects are defined and documented in the GenomicAlignments package.
See ?`inter-range-methods`
for
inter range transformations of a GenomicRanges or
GRangesList object.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 | ## S4 method for signature 'GenomicRanges'
shift(x, shift=0L, use.names=TRUE)
## S4 method for signature 'GenomicRanges'
narrow(x, start=NA, end=NA, width=NA, use.names=TRUE)
## S4 method for signature 'GenomicRanges'
resize(x, width, fix="start", use.names=TRUE, ignore.strand=FALSE)
## S4 method for signature 'GenomicRanges'
flank(x, width, start=TRUE, both=FALSE, use.names=TRUE,
ignore.strand=FALSE)
## S4 method for signature 'GenomicRanges'
promoters(x, upstream=2000, downstream=200, use.names=TRUE)
## S4 method for signature 'GenomicRanges'
restrict(x, start=NA, end=NA, keep.all.ranges=FALSE, use.names=TRUE)
## S4 method for signature 'GenomicRanges'
trim(x, use.names=TRUE)
|
x |
A GenomicRanges object. |
shift, use.names, start, end, width, both, fix, keep.all.ranges,
upstream, downstream |
See |
ignore.strand |
|
shift
behaves like the shift
method for
IntegerRanges objects. See
?`intra-range-methods`
for the details.
narrow
on a GenomicRanges object behaves
like on an IntegerRanges object. See
?`intra-range-methods`
for the details.
A major difference though is that it returns a GenomicRanges
object instead of an IntegerRanges object.
The returned object is parallel (i.e. same length and names)
to the original object x
.
resize
returns an object of the same type and length as
x
containing intervals that have been resized to width
width
based on the strand(x)
values. Elements where
strand(x) == "+"
or strand(x) == "*"
are anchored at
start(x)
and elements where strand(x) == "-"
are anchored
at the end(x)
. The use.names
argument determines whether
or not to keep the names on the ranges.
flank
returns an object of the same type and length
as x
containing intervals of width width
that flank
the intervals in x
. The start
argument takes a
logical indicating whether x
should be flanked at the
"start" (TRUE
) or the "end" (FALSE
), which for
strand(x) != "-"
is start(x)
and end(x)
respectively and for strand(x) == "-"
is end(x)
and
start(x)
respectively. The both
argument takes a
single logical value indicating whether the flanking region
width
positions extends into the range. If
both=TRUE
, the resulting range thus straddles the end
point, with width
positions on either side.
promoters
returns an object of the same type and length
as x
containing promoter ranges. Promoter ranges extend
around the transcription start site (TSS) which is defined as
start(x)
for ranges on the +
or *
strand
and as end(x)
for ranges on the -
strand.
The upstream
and downstream
arguments define the
number of nucleotides in the 5' and 3' direction, respectively.
More precisely, the output range is defined as
(start(x) - upstream) to (start(x) + downstream - 1)
for ranges on the +
or *
strand, and as
(end(x) - downstream + 1) to (end(x) + upstream)
for ranges on the -
strand.
Note that the returned object might contain out-of-bound ranges i.e. ranges that start before the first nucleotide position and/or end after the last nucleotide position of the underlying sequence.
restrict
returns an object of the same type and length as
x
containing restricted ranges for distinct seqnames. The
start
and end
arguments can be a named numeric vector of
seqnames for the ranges to be resticted or a numeric vector or length 1
if the restriction operation is to be applied to all the sequences in
x
. See ?`intra-range-methods`
for more
information about range restriction and for a description of the optional
arguments.
trim
trims out-of-bound ranges located on non-circular
sequences whose length is not NA.
P. Aboyoun and V. Obenchain
GenomicRanges, GRanges, and GRangesList objects.
The intra-range-methods man page in the IRanges package.
The IntegerRanges class in the IRanges package.
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 | ## ---------------------------------------------------------------------
## A. ON A GRanges OBJECT
## ---------------------------------------------------------------------
gr <- GRanges(
seqnames=Rle(paste("chr", c(1, 2, 1, 3), sep=""), c(1, 3, 2, 4)),
ranges=IRanges(1:10, width=10:1, names=letters[1:10]),
strand=Rle(strand(c("-", "+", "*", "+", "-")), c(1, 2, 2, 3, 2)),
score=1:10,
GC=seq(1, 0, length=10)
)
gr
shift(gr, 1)
narrow(gr[-10], start=2, end=-2)
resize(gr, width=10)
flank(gr, width=10)
restrict(gr, start=3, end=7)
gr <- GRanges("chr1", IRanges(rep(10, 3), width=6), c("+", "-", "*"))
promoters(gr, 2, 2)
## ---------------------------------------------------------------------
## B. ON A GRangesList OBJECT
## ---------------------------------------------------------------------
gr1 <- GRanges("chr2", IRanges(3, 6))
gr2 <- GRanges(c("chr1", "chr1"), IRanges(c(7,13), width=3),
strand=c("+", "-"))
gr3 <- GRanges(c("chr1", "chr2"), IRanges(c(1, 4), c(3, 9)),
strand="-")
grl <- GRangesList(gr1= gr1, gr2=gr2, gr3=gr3)
grl
resize(grl, width=20)
flank(grl, width=20)
restrict(grl, start=3)
|
Loading required package: stats4
Loading required package: BiocGenerics
Loading required package: parallel
Attaching package: ‘BiocGenerics’
The following objects are masked from ‘package:parallel’:
clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
clusterExport, clusterMap, parApply, parCapply, parLapply,
parLapplyLB, parRapply, parSapply, parSapplyLB
The following objects are masked from ‘package:stats’:
IQR, mad, sd, var, xtabs
The following objects are masked from ‘package:base’:
anyDuplicated, append, as.data.frame, basename, cbind, colnames,
dirname, do.call, duplicated, eval, evalq, Filter, Find, get, grep,
grepl, intersect, is.unsorted, lapply, Map, mapply, match, mget,
order, paste, pmax, pmax.int, pmin, pmin.int, Position, rank,
rbind, Reduce, rownames, sapply, setdiff, sort, table, tapply,
union, unique, unsplit, which.max, which.min
Loading required package: S4Vectors
Attaching package: ‘S4Vectors’
The following object is masked from ‘package:base’:
expand.grid
Loading required package: IRanges
Loading required package: GenomeInfoDb
GRanges object with 10 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
a chr1 1-10 - | 1 1.000000
b chr2 2-10 + | 2 0.888889
c chr2 3-10 + | 3 0.777778
d chr2 4-10 * | 4 0.666667
e chr1 5-10 * | 5 0.555556
f chr1 6-10 + | 6 0.444444
g chr3 7-10 + | 7 0.333333
h chr3 8-10 + | 8 0.222222
i chr3 9-10 - | 9 0.111111
j chr3 10 - | 10 0.000000
-------
seqinfo: 3 sequences from an unspecified genome; no seqlengths
GRanges object with 10 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
a chr1 2-11 - | 1 1.000000
b chr2 3-11 + | 2 0.888889
c chr2 4-11 + | 3 0.777778
d chr2 5-11 * | 4 0.666667
e chr1 6-11 * | 5 0.555556
f chr1 7-11 + | 6 0.444444
g chr3 8-11 + | 7 0.333333
h chr3 9-11 + | 8 0.222222
i chr3 10-11 - | 9 0.111111
j chr3 11 - | 10 0.000000
-------
seqinfo: 3 sequences from an unspecified genome; no seqlengths
GRanges object with 9 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
a chr1 2-9 - | 1 1.000000
b chr2 3-9 + | 2 0.888889
c chr2 4-9 + | 3 0.777778
d chr2 5-9 * | 4 0.666667
e chr1 6-9 * | 5 0.555556
f chr1 7-9 + | 6 0.444444
g chr3 8-9 + | 7 0.333333
h chr3 9 + | 8 0.222222
i chr3 10-9 - | 9 0.111111
-------
seqinfo: 3 sequences from an unspecified genome; no seqlengths
GRanges object with 10 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
a chr1 1-10 - | 1 1.000000
b chr2 2-11 + | 2 0.888889
c chr2 3-12 + | 3 0.777778
d chr2 4-13 * | 4 0.666667
e chr1 5-14 * | 5 0.555556
f chr1 6-15 + | 6 0.444444
g chr3 7-16 + | 7 0.333333
h chr3 8-17 + | 8 0.222222
i chr3 1-10 - | 9 0.111111
j chr3 1-10 - | 10 0.000000
-------
seqinfo: 3 sequences from an unspecified genome; no seqlengths
GRanges object with 10 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
a chr1 11-20 - | 1 1.000000
b chr2 -8-1 + | 2 0.888889
c chr2 -7-2 + | 3 0.777778
d chr2 -6-3 * | 4 0.666667
e chr1 -5-4 * | 5 0.555556
f chr1 -4-5 + | 6 0.444444
g chr3 -3-6 + | 7 0.333333
h chr3 -2-7 + | 8 0.222222
i chr3 11-20 - | 9 0.111111
j chr3 11-20 - | 10 0.000000
-------
seqinfo: 3 sequences from an unspecified genome; no seqlengths
GRanges object with 8 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
a chr1 3-7 - | 1 1.000000
b chr2 3-7 + | 2 0.888889
c chr2 3-7 + | 3 0.777778
d chr2 4-7 * | 4 0.666667
e chr1 5-7 * | 5 0.555556
f chr1 6-7 + | 6 0.444444
g chr3 7 + | 7 0.333333
h chr3 8-7 + | 8 0.222222
-------
seqinfo: 3 sequences from an unspecified genome; no seqlengths
GRanges object with 3 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 8-11 +
[2] chr1 14-17 -
[3] chr1 8-11 *
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr2 3-6 *
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr2
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 7-9 +
[2] chr1 13-15 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr3
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 1-3 -
[2] chr2 4-9 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr2 3-22 *
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr2
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 7-26 +
[2] chr1 -4-15 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr3
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 -16-3 -
[2] chr2 -10-9 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr2 -17-2 *
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr2
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 -13-6 +
[2] chr1 16-35 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr3
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 4-23 -
[2] chr2 10-29 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr2 3-6 *
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr2
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 7-9 +
[2] chr1 13-15 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
$gr3
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr1 3 -
[2] chr2 4-9 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.