BSgenomeViews-class {BSgenome} | R Documentation |
BSgenomeViews objects
Description
The BSgenomeViews class is a container for storing a set of genomic positions on a BSgenome object, called the "subject" in this context.
Usage
## Constructor
## ------------
BSgenomeViews(subject, granges)
## Accessors
## ---------
## S4 method for signature 'BSgenomeViews'
subject(x)
## S4 method for signature 'BSgenomeViews'
granges(x, use.mcols=FALSE)
## S4 method for signature 'BSgenomeViews'
length(x)
## S4 method for signature 'BSgenomeViews'
names(x)
## S4 method for signature 'BSgenomeViews'
seqnames(x)
## S4 method for signature 'BSgenomeViews'
start(x)
## S4 method for signature 'BSgenomeViews'
end(x)
## S4 method for signature 'BSgenomeViews'
width(x)
## S4 method for signature 'BSgenomeViews'
strand(x)
## S4 method for signature 'BSgenomeViews'
ranges(x, use.mcols=FALSE)
## S4 method for signature 'BSgenomeViews'
elementNROWS(x)
## S4 method for signature 'BSgenomeViews'
seqinfo(x)
## DNAStringSet methods
## --------------------
## S4 method for signature 'BSgenomeViews'
seqtype(x)
## S4 method for signature 'BSgenomeViews'
nchar(x, type="chars", allowNA=FALSE)
## S4 method for signature 'BSgenomeViews'
unlist(x, recursive=TRUE, use.names=TRUE)
## S4 method for signature 'BSgenomeViews'
alphabetFrequency(x, as.prob=FALSE, collapse=FALSE, baseOnly=FALSE)
## S4 method for signature 'BSgenomeViews'
hasOnlyBaseLetters(x)
## S4 method for signature 'BSgenomeViews'
uniqueLetters(x)
## S4 method for signature 'BSgenomeViews'
letterFrequency(x, letters, OR="|", as.prob=FALSE, collapse=FALSE)
## S4 method for signature 'BSgenomeViews'
oligonucleotideFrequency(x, width, step=1,
as.prob=FALSE, as.array=FALSE,
fast.moving.side="right", with.labels=TRUE, simplify.as="matrix")
## S4 method for signature 'BSgenomeViews'
nucleotideFrequencyAt(x, at, as.prob=FALSE, as.array=TRUE,
fast.moving.side="right", with.labels=TRUE)
## S4 method for signature 'BSgenomeViews'
consensusMatrix(x, as.prob=FALSE, shift=0L, width=NULL, baseOnly=FALSE)
## S4 method for signature 'BSgenomeViews'
consensusString(x, ambiguityMap=IUPAC_CODE_MAP, threshold=0.25,
shift=0L, width=NULL)
Arguments
subject |
A BSgenome object or the name of a reference genome specified
in a way that is accepted by the |
granges |
A GRanges object containing ranges relative to
the genomic sequences stored in |
x |
A BSgenomeViews object. |
use.mcols |
|
type , allowNA , recursive , use.names |
Ignored. |
as.prob , letters , OR , width |
See |
collapse , baseOnly |
See |
step , as.array , fast.moving.side , with.labels , simplify.as , at |
See |
shift , ambiguityMap , threshold |
See |
Constructors
-
BSgenomeViews(subject, granges)
: Make a BSgenomeViews object by putting the views specified bygranges
on top of the genomic sequences stored insubject
. See above for how argumentsubject
andgranges
should be specified. -
Views(subject, granges)
: Equivalent toBSgenomeViews(subject, granges)
. Provided for convenience.
Accessors
In the code snippets below, x
is a BSgenomeViews object.
-
subject(x)
: Return the BSgenome object containing the full genomic sequences on top of which the views inx
are defined. -
granges(x, use.mcols=FALSE)
: Return the genomic ranges of the views as a GRanges object. These ranges are relative to the genomic sequences stored insubject(x)
. -
length(x)
: The number of views inx
. -
names(x)
: The names of the views inx
. -
seqnames(x)
,start(x)
,end(x)
,width(x)
,strand(x)
: Equivalent toseqnames(granges(x))
,start(granges(x))
,end(granges(x))
,width(granges(x))
,strand(granges(x))
, respectively. -
ranges(x, use.mcols=FALSE)
: Equivalent toranges(granges(x, use.mcols), use.mcols)
. -
elementNROWS(x)
: Equivalent towidth(x)
. -
seqinfo(x)
: Equivalent toseqinfo(subject(x))
and toseqinfo(granges(x))
(both are guaranteed to be the same). See?seqinfo
in the GenomeInfoDb package for more information.
Coercion
In the code snippets below, x
is a BSgenomeViews object.
-
as(x, "DNAStringSet")
: Turnx
into a DNAStringSet object by extxracting the DNA sequence corresponding to each view. Alternativelyas(x, "XStringSet")
can be used for this, and is equivalent toas(x, "DNAStringSet")
. -
as.character(x)
: Equivalent toas.character(as(x, "DNAStringSet"))
. -
as.data.frame(x)
: Turnx
into a data.frame.
Subsetting
-
x[i]
: Select the views specified byi
. -
x[[i]]
: Extract the one view specified byi
.
DNAStringSet methods
For convenience, some methods defined for DNAStringSet
objects in the Biostrings package can be used directly on a
BSgenomeViews object. In that case, everything happens like if the
BSgenomeViews object x
was turned into a
DNAStringSet object (with as(x, "DNAStringSet")
)
before it's passed to the method for DNAStringSet objects.
At the moment, the list of such methods is:
seqtype
,
nchar,XStringSet-method
,
unlist,XStringSet-method
,
alphabetFrequency
,
hasOnlyBaseLetters
,
uniqueLetters
,
letterFrequency
,
oligonucleotideFrequency
,
nucleotideFrequencyAt
,
consensusMatrix
,
and consensusString
.
See the corresponding man page in the Biostrings package for a description of these methods.
Author(s)
H. Pagès
See Also
The BSgenome class.
The GRanges class in the GenomicRanges package.
The DNAStringSet class in the Biostrings package.
The seqinfo and related getters in the GenomeInfoDb package for getting the sequence information stored in an object.
-
TxDb objects in the GenomicFeatures package.
Examples
library(BSgenome.Mmusculus.UCSC.mm10)
genome <- BSgenome.Mmusculus.UCSC.mm10
library(TxDb.Mmusculus.UCSC.mm10.knownGene)
txdb <- TxDb.Mmusculus.UCSC.mm10.knownGene
ex <- exons(txdb, columns=c("exon_id", "tx_name", "gene_id"))
v <- Views(genome, ex)
v
subject(v)
granges(v)
seqinfo(v)
as(v, "DNAStringSet")
v10 <- v[1:10] # select the first 10 views
subject(v10) # same as subject(v)
granges(v10)
seqinfo(v10) # same as seqinfo(v)
as(v10, "DNAStringSet")
alphabetFrequency(v10)
alphabetFrequency(v10, collapse=TRUE)
v12 <- v[width(v) <= 12] # select the views of 12 nucleotides or less
head(as.data.frame(v12))
trinucleotideFrequency(v12, simplify.as="collapsed")
## BSgenomeViews objects are list-like objects. That is, the
## BSgenomeViews class derives from List and typical list/List
## operations (e.g. [[, elementNROWS(), unlist(), elementType(),
## etc...) work on these objects:
is(v12, "List") # TRUE
v12[[2]]
head(elementNROWS(v12)) # elementNROWS(v) is the same as width(v)
unlist(v12)
elementType(v12)