The Hits class is a container for representing a set of hits between
a set of left nodes and a set of right nodes.
Note that only the hits are stored in the object. No information about
the left or right nodes is stored, except their number.
For example, the findOverlaps function, defined
and documented in the IRanges package, returns the hits between
the query and subject arguments in a Hits object.
2 integer vectors of the same length.
The values in from must be >= 1 and <= nLnode.
The values in to must be >= 1 and <= nRnode.
nLnode, nRnode
Number of left and right nodes.
...
Metadata columns to set on the Hits object. All the metadata columns must
be vector-like objects of the same length as from and to.
sort.by.query
Should the hits in the returned object be sorted by query? If yes, then
a SortedByQueryHits object is returned (SortedByQueryHits is a subclass
of Hits).
nnode
Number of nodes.
Accessors
In the code snippets below, x is a Hits object.
length(x): get the number of hits
from(x): Equivalent to as.data.frame(x)[[1]].
to(x): Equivalent to as.data.frame(x)[[2]].
nLnode(x), nrow(x): get the number of left nodes
nRnode(x), ncol(x): get the number of right nodes
countLnodeHits(x): Counts the number of hits for
each left node, returning an integer vector.
countRnodeHits(x): Counts the number of hits for
each right node, returning an integer vector.
The following accessors are just aliases for the above accessors:
queryHits(x): alias for from(x).
subjectHits(x): alias for to(x).
queryLength(x): alias for nLnode(x).
subjectLength(x): alias for nRnode(x).
countQueryHits(x): alias for countLnodeHits(x).
countSubjectHits(x): alias for countRnodeHits(x).
Coercion
In the code snippets below, x is a Hits object.
as.matrix(x): Coerces x to a two
column integer matrix, with each row representing a hit
between a left node (first column) and a right node (second
column).
as.table(x): Counts the number of hits for
each left node in x and outputs the counts as a table.
Subsetting
In the code snippets below, x is a Hits object.
x[i]: Subset the Hits object.
Other transformations
In the code snippets below, x is a Hits object.
t(x):
Transpose x by interchanging the left and right nodes. This
allows, for example, counting the number of hits for each right node
using as.table.
remapHits(x, Lnodes.remapping=NULL, new.nLnode=NA,
Rnodes.remapping=NULL, new.nRnode=NA):
Only supports SortedByQueryHits objects at the moment.
Remaps the left and/or right nodes in x. The left nodes are
remapped thru the map specified via the Lnodes.remapping and
new.nLnode arguments. The right nodes are remapped thru the
map specified via the Rnodes.remapping and new.nRnode
arguments.
Lnodes.remapping must represent a function defined on the
1..M interval that takes values in the 1..N interval, where N is
nLnode(x) and M is the value specified by the user via the
new.nLnode argument. Note that this mapping function doesn't
need to be injective or surjective. Also it is not represented by an R
function but by an integer vector of length M with no NAs. More precisely
Lnodes.remapping can be NULL (identity map), or a vector of
nLnode(x) non-NA integers that are >= 1 and
<= new.nLnode, or a factor of length nLnode(x)
with no NAs (a factor is treated as an integer vector, and, if missing,
new.nLnode is taken to be its number of levels). Note that
a factor will typically be used to represent a mapping function that is
not injective.
The same applies to the Rnodes.remapping.
remapHits returns a Hits object where from(x) and
to(x) have been remapped thru the 2 specified maps. This
remapping is actually only the 1st step of the transformation, and is
followed by 2 additional steps: (2) the removal of duplicated hits,
and (3) the reordering of the hits (first by query hits, then by subject
hits). Note that if the 2 maps are injective then the remapping won't
introduce duplicated hits, so, in that case, step (2) is a no-op (but
is still performed). Also if the "query map" is strictly ascending and
the "subject map" ascending then the remapping will preserve the order
of the hits, so, in that case, step (3) is also a no-op (but is still
performed).
SelfHits
A SelfHits object is a Hits object where the left and right nodes are
identical. For a SelfHits object x, nLnode(x) is equal to
nRnode(x). The object can be seen as an oriented graph where
nLnode is the nb of nodes and the hits are the (oriented) edges.
SelfHits objects support the same set of accessors as Hits objects
plus the nnode() accessor that is equivalent to nLnode()
and nRnode().
We also provide two little utilities to operate on a SelfHits object
x:
isSelfHit(x): A self hit is an edge from a node to
itself. isSelfHit(x) returns a logical vector parallel
to x indicating which elements in x are self hits.
isRedundantHit(x): When there is more than 1 edge between
2 given nodes (regardless of orientation), the extra edges are considered
to be redundant hits. isRedundantHit(x) returns a logical
vector parallel to x indicating which elements in x
are redundant hits.
Author(s)
Michael Lawrence and Herv<c3><83><c2><a9> Pag<c3><83><c2><a8>s
See Also
Hits-comparison for comparing and ordering hits.
The findOverlaps function in the
IRanges package which returns SortedByQueryHits object.
Hits-examples in the IRanges package, for
some examples of Hits object basic manipulation.
setops-methods in the IRanges package,
for set operations on Hits objects.
Examples
from <- c(5, 2, 3, 3, 3, 2)
to <- c(11, 15, 5, 4, 5, 11)
id <- letters[1:6]
Hits(from, to, 7, 15, id)
Hits(from, to, 7, 15, id, sort.by.query=TRUE)
## ---------------------------------------------------------------------
## selectHits()
## ---------------------------------------------------------------------
x <- c("a", "b", "a", "c", "d")
table <- c("a", "e", "d", "a", "a", "d")
hits <- findMatches(x, table) # sorts the hits by query
hits
selectHits(hits, select="all") # no-op
selectHits(hits, select="first")
selectHits(hits, select="last")
selectHits(hits, select="arbitrary")
selectHits(hits, select="count")
## ---------------------------------------------------------------------
## remapHits()
## ---------------------------------------------------------------------
Lnodes.remapping <- factor(c(a="A", b="B", c="C", d="D")[x],
levels=LETTERS[1:4])
remapHits(hits, Lnodes.remapping=Lnodes.remapping)
## See ?`Hits-examples` in the IRanges package for more examples of basic
## manipulation of Hits objects.
## ---------------------------------------------------------------------
## SelfHits objects
## ---------------------------------------------------------------------
hits2 <- SelfHits(c(2, 3, 3, 3, 3, 3, 4, 4, 4), c(4, 3, 2:4, 2, 2:3, 2), 4)
## Hits 2 and 4 are self hits (from 3rd node to itself):
which(isSelfHit(hits2))
## Hits 4, 6, 7, 8, and 9, are redundant hits:
which(isRedundantHit(hits2))
hits3 <- findMatches(x)
hits3[!isSelfHit(hits3)]
hits3[!(isSelfHit(hits3) | isRedundantHit(hits3))]
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> library(S4Vectors)
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, xtabs
The following objects are masked from 'package:base':
Filter, Find, Map, Position, Reduce, anyDuplicated, append,
as.data.frame, cbind, colnames, 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, rownames, sapply, setdiff, sort, table, tapply, union,
unique, unsplit
Attaching package: 'S4Vectors'
The following objects are masked from 'package:base':
colMeans, colSums, expand.grid, rowMeans, rowSums
> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/S4Vectors/Hits-class.Rd_%03d_medium.png", width=480, height=480)
> ### Name: Hits-class
> ### Title: Hits objects
> ### Aliases: class:Hits Hits-class Hits class:SelfHits SelfHits-class
> ### SelfHits class:SortedByQueryHits SortedByQueryHits-class
> ### SortedByQueryHits class:SortedByQuerySelfHits
> ### SortedByQuerySelfHits-class SortedByQuerySelfHits
> ### parallelSlotNames,Hits-method from from,Hits-method to to,Hits-method
> ### nLnode nLnode,Hits-method nRnode nRnode,Hits-method nnode
> ### nnode,SelfHits-method countLnodeHits countLnodeHits,Hits-method
> ### countRnodeHits countRnodeHits,Hits-method queryHits subjectHits
> ### queryLength subjectLength countQueryHits countSubjectHits
> ### updateObject,Hits-method coerce,Hits,SortedByQueryHits-method
> ### as.matrix,Hits-method as.table,Hits-method
> ### extractROWS,classNameForDisplay,ANY-method
> ### classNameForDisplay,Hits-method classNameForDisplay,SelfHits-method
> ### show,Hits-method c,Hits-method t,Hits-method remapHits isSelfHit
> ### isRedundantHit
> ### Keywords: methods classes
>
> ### ** Examples
>
> from <- c(5, 2, 3, 3, 3, 2)
> to <- c(11, 15, 5, 4, 5, 11)
> id <- letters[1:6]
>
> Hits(from, to, 7, 15, id)
Hits object with 6 hits and 1 metadata column:
from to | id
<integer> <integer> | <character>
[1] 5 11 | a
[2] 2 15 | b
[3] 3 5 | c
[4] 3 4 | d
[5] 3 5 | e
[6] 2 11 | f
-------
nLnode: 7 / nRnode: 15
> Hits(from, to, 7, 15, id, sort.by.query=TRUE)
Hits object with 6 hits and 1 metadata column:
queryHits subjectHits | id
<integer> <integer> | <character>
[1] 2 15 | b
[2] 2 11 | f
[3] 3 5 | c
[4] 3 4 | d
[5] 3 5 | e
[6] 5 11 | a
-------
queryLength: 7 / subjectLength: 15
>
> ## ---------------------------------------------------------------------
> ## selectHits()
> ## ---------------------------------------------------------------------
>
> x <- c("a", "b", "a", "c", "d")
> table <- c("a", "e", "d", "a", "a", "d")
> hits <- findMatches(x, table) # sorts the hits by query
> hits
Hits object with 8 hits and 0 metadata columns:
queryHits subjectHits
<integer> <integer>
[1] 1 1
[2] 1 4
[3] 1 5
[4] 3 1
[5] 3 4
[6] 3 5
[7] 5 3
[8] 5 6
-------
queryLength: 5 / subjectLength: 6
>
> selectHits(hits, select="all") # no-op
Hits object with 8 hits and 0 metadata columns:
queryHits subjectHits
<integer> <integer>
[1] 1 1
[2] 1 4
[3] 1 5
[4] 3 1
[5] 3 4
[6] 3 5
[7] 5 3
[8] 5 6
-------
queryLength: 5 / subjectLength: 6
> selectHits(hits, select="first")
[1] 1 NA 1 NA 3
> selectHits(hits, select="last")
[1] 5 NA 5 NA 6
> selectHits(hits, select="arbitrary")
[1] 5 NA 5 NA 6
> selectHits(hits, select="count")
[1] 3 0 3 0 2
>
> ## ---------------------------------------------------------------------
> ## remapHits()
> ## ---------------------------------------------------------------------
>
> Lnodes.remapping <- factor(c(a="A", b="B", c="C", d="D")[x],
+ levels=LETTERS[1:4])
> remapHits(hits, Lnodes.remapping=Lnodes.remapping)
Hits object with 5 hits and 0 metadata columns:
queryHits subjectHits
<integer> <integer>
[1] 1 1
[2] 1 4
[3] 1 5
[4] 4 3
[5] 4 6
-------
queryLength: 4 / subjectLength: 6
>
> ## See ?`Hits-examples` in the IRanges package for more examples of basic
> ## manipulation of Hits objects.
>
> ## ---------------------------------------------------------------------
> ## SelfHits objects
> ## ---------------------------------------------------------------------
>
> hits2 <- SelfHits(c(2, 3, 3, 3, 3, 3, 4, 4, 4), c(4, 3, 2:4, 2, 2:3, 2), 4)
> ## Hits 2 and 4 are self hits (from 3rd node to itself):
> which(isSelfHit(hits2))
[1] 2 4
> ## Hits 4, 6, 7, 8, and 9, are redundant hits:
> which(isRedundantHit(hits2))
[1] 4 6 7 8 9
>
> hits3 <- findMatches(x)
> hits3[!isSelfHit(hits3)]
SelfHits object with 2 hits and 0 metadata columns:
queryHits subjectHits
<integer> <integer>
[1] 1 3
[2] 3 1
-------
queryLength: 5 / subjectLength: 5
> hits3[!(isSelfHit(hits3) | isRedundantHit(hits3))]
SelfHits object with 1 hit and 0 metadata columns:
queryHits subjectHits
<integer> <integer>
[1] 1 3
-------
queryLength: 5 / subjectLength: 5
>
>
>
>
>
> dev.off()
null device
1
>