Description Constructors Accessors Coercion Subsetting Combining Looping Author(s) See Also Examples
The GRangesList class is a container for storing a collection of GRanges objects. It is a subclass of GenomicRangesList. It exists in 2 flavors: SimpleGRangesList and CompressedGRangesList. Each flavor uses a particular internal representation. The CompressedGRangesList flavor is the default. It is particularly efficient for storing a large number of list elements and operating on them.
GRangesList(..., compress=TRUE)
:Creates a GRangesList object using the GRanges objects
supplied in ...
, either consecutively or in a list.
By default a CompressedGRangesList instance is returned,
that is, a GRangesList object of the CompressedGRangesList flavor.
Use compress=FALSE
to get a SimpleGRangesList instance
instead.
makeGRangesListFromFeatureFragments(seqnames=Rle(factor()),
fragmentStarts=list(), fragmentEnds=list(),
fragmentWidths=list(),
strand=character(0), sep=",")
:Constructs a GRangesList object from a list of fragmented features. See the Examples section below.
In the code snippets below, x
is a GRangesList object.
length(x)
:Get the number of list elements.
names(x)
, names(x) <- value
:Get or set the names on x
.
seqnames(x)
, seqnames(x) <- value
:Get or set the sequence names in the form of an RleList.
value
can be an RleList or CharacterList object.
ranges(x, use.mcols=FALSE)
, ranges(x) <- value
:Get or set the ranges in the form of a CompressedIRangesList.
value
can be an IntegerRangesList object.
start(x)
, start(x) <- value
:Get or set start(ranges(x))
.
end(x)
, end(x) <- value
:Get or set end(ranges(x))
.
width(x)
, width(x) <- value
:Get or set width(ranges(x))
.
strand(x)
, strand(x) <- value
:Get or set the strand in the form of an RleList object.
value
can be an RleList, CharacterList or
single character. value
as a single character converts all
ranges in x
to the same value
; for selective strand
conversion (i.e., mixed +
and -
) use RleList
or CharacterList.
mcols(x, use.names=FALSE)
, mcols(x) <- value
:Get or set the metadata columns.
value
can be NULL
, or a data.frame-like object (i.e.
DataFrame or data.frame) holding element-wise metadata.
elementNROWS(x)
:Get a vector of the length
of each of the list element.
isEmpty(x)
:Returns a logical indicating either if the GRangesList has no elements or if all its elements are empty.
seqinfo(x)
, seqinfo(x) <- value
:Get or set the information about the underlying sequences.
value
must be a Seqinfo object.
seqlevels(x)
,
seqlevels(x, pruning.mode=c("error", "coarse", "fine", "tidy"))
<- value
:Get or set the sequence levels.
seqlevels(x)
is equivalent to seqlevels(seqinfo(x))
or to levels(seqnames(x))
, those 2 expressions being
guaranteed to return identical character vectors on a GRangesList
object. value
must be a character vector with no NAs.
See ?seqlevels
for more information.
seqlengths(x)
, seqlengths(x) <- value
:Get or set the sequence lengths.
seqlengths(x)
is equivalent to seqlengths(seqinfo(x))
.
value
can be a named non-negative integer or numeric vector
eventually with NAs.
isCircular(x)
, isCircular(x) <- value
:Get or set the circularity flags.
isCircular(x)
is equivalent to isCircular(seqinfo(x))
.
value
must be a named logical vector eventually with NAs.
genome(x)
, genome(x) <- value
:Get or set the genome identifier or assembly name for each sequence.
genome(x)
is equivalent to genome(seqinfo(x))
.
value
must be a named character vector eventually with NAs.
seqlevelsStyle(x)
, seqlevelsStyle(x) <- value
:Get or set the seqname style for x
.
See the seqlevelsStyle generic getter and setter
in the GenomeInfoDb package for more information.
score(x), score(x) <- value
:Get or set the score
metadata column.
In the code snippets below, x
is a GRangesList object.
as.data.frame(x, row.names=NULL, optional=FALSE, ...,
value.name="value", use.outer.mcols=FALSE,
group_name.as.factor=FALSE)
:Coerces x
to a data.frame
. See as.data.frame on the
List
man page for details (?List
).
as.list(x, use.names = TRUE)
:Creates a list containing the elements of x
.
as(x, "IRangesList")
:Turns x
into an IRangesList object.
When x
is a list
of GRanges
, it can be coerced to a
GRangesList
.
as(x, "GRangesList")
:Turns x
into a GRangesList
.
In the following code snippets, x
is a GRangesList object.
x[i, j]
, x[i, j] <- value
:Get or set elements i
with optional metadata columns
mcols(x)[,j]
, where i
can be missing; an NA-free
logical, numeric, or character vector; a logical-Rle object,
or an AtomicList object.
x[[i]]
, x[[i]] <- value
:Get or set element i
, where i
is a numeric or character
vector of length 1.
x$name
, x$name <- value
:Get or set element name
, where name
is a name or character
vector of length 1.
head(x, n = 6L)
:If n
is non-negative, returns the first n elements of the
GRangesList object.
If n
is negative, returns all but the last abs(n)
elements
of the GRangesList object.
rep(x, times, length.out, each)
:Repeats the values in x
through one of the following conventions:
times
Vector giving the number of times to repeat each
element if of length length(x)
, or to repeat the whole vector
if of length 1.
length.out
Non-negative integer. The desired length of the output vector.
each
Non-negative integer. Each element of x
is
repeated each
times.
subset(x, subset)
:Returns a new object of the same class as x
made of the subset
using logical vector subset
, where missing values are taken as
FALSE
.
tail(x, n = 6L)
:If n
is non-negative, returns the last n list elements of the
GRangesList object.
If n
is negative, returns all but the first abs(n)
list elements of the GRangesList object.
In the code snippets below, x
is a GRangesList object.
c(x, ...)
:Combines x
and the GRangesList objects in ...
together. Any object in ...
must belong to the same class
as x
, or to one of its subclasses, or must be NULL
.
The result is an object of the same class as x
.
append(x, values, after = length(x))
:Inserts the values
into x
at the position given by
after
, where x
and values
are of the same
class.
unlist(x, recursive = TRUE, use.names = TRUE)
:Concatenates the elements of x
into a single GRanges
object.
In the code snippets below, x
is a GRangesList object.
endoapply(X, FUN, ...)
:Similar to lapply
, but performs an endomorphism,
i.e. returns an object of class(X)
.
lapply(X, FUN, ...)
:Like the standard lapply
function defined in the
base package, the lapply
method for GRangesList objects
returns a list of the same length as X
, with each element being
the result of applying FUN
to the corresponding element of
X
.
Map(f, ...)
:Applies a function to the corresponding elements of given GRangesList objects.
mapply(FUN, ..., MoreArgs=NULL, SIMPLIFY=TRUE, USE.NAMES=TRUE)
:Like the standard mapply
function defined in the
base package, the mapply
method for GRangesList objects is a
multivariate version of sapply
.
mendoapply(FUN, ..., MoreArgs = NULL)
:Similar to mapply
, but performs an endomorphism
across multiple objects, i.e. returns an object of
class(list(...)[[1]])
.
Reduce(f, x, init, right = FALSE, accumulate = FALSE)
:Uses a binary function to successively combine the elements of x
and a possibly given initial value.
f
A binary argument function.
init
An R object of the same kind as the elements of x
.
right
A logical indicating whether to proceed from left to right (default) or from right to left.
nomatch
The value to be returned in the case when "no match" (no element satisfying the predicate) is found.
sapply(X, FUN, ..., simplify=TRUE, USE.NAMES=TRUE)
:Like the standard sapply
function defined in
the base package, the sapply
method for GRangesList objects
is a user-friendly version of lapply
by default returning a vector
or matrix if appropriate.
P. Aboyoun & H. Pag<c3><a8>s
GRanges objects.
seqinfo
in the GenomeInfoDb
package.
IntegerRangesList objects in the IRanges package.
RleList objects in the IRanges package.
DataFrameList objects in the IRanges package.
intra-range-methods, inter-range-methods, coverage-methods, setops-methods, and findOverlaps-methods.
GenomicRangesList objects.
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 | ## Construction with GRangesList():
gr1 <- GRanges("chr2", IRanges(3, 6),
strand="+", score=5L, GC=0.45)
gr2 <- GRanges(c("chr1", "chr1"), IRanges(c(7,13), width=3),
strand=c("+", "-"), score=3:4, GC=c(0.3, 0.5))
gr3 <- GRanges(c("chr1", "chr2"), IRanges(c(1, 4), c(3, 9)),
strand=c("-", "-"), score=c(6L, 2L), GC=c(0.4, 0.1))
grl <- GRangesList(gr1=gr1, gr2=gr2, gr3=gr3)
grl
## Summarizing elements:
elementNROWS(grl)
table(seqnames(grl))
## Extracting subsets:
grl[seqnames(grl) == "chr1", ]
grl[seqnames(grl) == "chr1" & strand(grl) == "+", ]
## Renaming the underlying sequences:
seqlevels(grl)
seqlevels(grl) <- sub("chr", "Chrom", seqlevels(grl))
grl
## Coerce to IRangesList (seqnames and strand information is lost):
as(grl, "IRangesList")
## isDisjoint():
isDisjoint(grl)
## disjoin():
disjoin(grl) # metadata columns and order NOT preserved
## Construction with makeGRangesListFromFeatureFragments():
filepath <- system.file("extdata", "feature_frags.txt",
package="GenomicRanges")
featfrags <- read.table(filepath, header=TRUE, stringsAsFactors=FALSE)
grl2 <- with(featfrags,
makeGRangesListFromFeatureFragments(seqnames=targetName,
fragmentStarts=targetStart,
fragmentWidths=blockSizes,
strand=strand))
names(grl2) <- featfrags$RefSeqID
grl2
|
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':
Filter, Find, Map, Position, Reduce, anyDuplicated, append,
as.data.frame, basename, cbind, colMeans, colSums, colnames,
dirname, do.call, duplicated, eval, evalq, get, grep, grepl,
intersect, is.unsorted, lapply, lengths, mapply, match, mget,
order, paste, pmax, pmax.int, pmin, pmin.int, rank, rbind,
rowMeans, rowSums, rownames, sapply, setdiff, sort, table, tapply,
union, unique, unsplit, which, 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
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
[1] chr2 3-6 + | 5 0.45
$gr2
GRanges object with 2 ranges and 2 metadata columns:
seqnames ranges strand | score GC
[1] chr1 7-9 + | 3 0.3
[2] chr1 13-15 - | 4 0.5
$gr3
GRanges object with 2 ranges and 2 metadata columns:
seqnames ranges strand | score GC
[1] chr1 1-3 - | 6 0.4
[2] chr2 4-9 - | 2 0.1
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
gr1 gr2 gr3
1 2 2
chr2 chr1
gr1 1 0
gr2 0 2
gr3 1 1
GRangesList object of length 3:
$gr1
GRanges object with 0 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
$gr2
GRanges object with 2 ranges and 2 metadata columns:
seqnames ranges strand | score GC
[1] chr1 7-9 + | 3 0.3
[2] chr1 13-15 - | 4 0.5
$gr3
GRanges object with 1 range and 2 metadata columns:
seqnames ranges strand | score GC
[1] chr1 1-3 - | 6 0.4
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
GRangesList object of length 3:
$gr1
GRanges object with 0 ranges and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
$gr2
GRanges object with 1 range and 2 metadata columns:
seqnames ranges strand | score GC
[1] chr1 7-9 + | 3 0.3
$gr3
GRanges object with 0 ranges and 2 metadata columns:
seqnames ranges strand | score GC
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
[1] "chr2" "chr1"
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 2 metadata columns:
seqnames ranges strand | score GC
<Rle> <IRanges> <Rle> | <integer> <numeric>
[1] Chrom2 3-6 + | 5 0.45
$gr2
GRanges object with 2 ranges and 2 metadata columns:
seqnames ranges strand | score GC
[1] Chrom1 7-9 + | 3 0.3
[2] Chrom1 13-15 - | 4 0.5
$gr3
GRanges object with 2 ranges and 2 metadata columns:
seqnames ranges strand | score GC
[1] Chrom1 1-3 - | 6 0.4
[2] Chrom2 4-9 - | 2 0.1
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
IRangesList of length 3
$gr1
IRanges object with 1 range and 2 metadata columns:
start end width | score GC
<integer> <integer> <integer> | <integer> <numeric>
[1] 3 6 4 | 5 0.45
$gr2
IRanges object with 2 ranges and 2 metadata columns:
start end width | score GC
<integer> <integer> <integer> | <integer> <numeric>
[1] 7 9 3 | 3 0.3
[2] 13 15 3 | 4 0.5
$gr3
IRanges object with 2 ranges and 2 metadata columns:
start end width | score GC
<integer> <integer> <integer> | <integer> <numeric>
[1] 1 3 3 | 6 0.4
[2] 4 9 6 | 2 0.1
gr1 gr2 gr3
TRUE TRUE TRUE
GRangesList object of length 3:
$gr1
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] Chrom2 3-6 +
$gr2
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
[1] Chrom1 7-9 +
[2] Chrom1 13-15 -
$gr3
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
[1] Chrom2 4-9 -
[2] Chrom1 1-3 -
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
GRangesList object of length 3:
$XM_001065892.1
GRanges object with 2 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chr4 124513961-124514087 +
[2] chr4 124514706-124515691 +
$XM_578205.2
GRanges object with 3 ranges and 0 metadata columns:
seqnames ranges strand
[1] chr2 155875533-155876067 -
[2] chr2 155879894-155880030 -
[3] chr2 155895543-155895690 -
$NM_012543.2
GRanges object with 4 ranges and 0 metadata columns:
seqnames ranges strand
[1] chr1 96173572-96174077 +
[2] chr1 96174920-96175330 +
[3] chr1 96176574-96176785 +
[4] chr1 96177991-96178484 +
-------
seqinfo: 3 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.