Create SummarizedExperiment representations by transforming ragged assays to rectangular form.
Source:R/coerce-functions.R
coerce-functions.RdThese methods transform RaggedExperiment
objects to similar SummarizedExperiment objects. They do
so by transforming assay data to more rectangular
representations, following the rules outlined for similarly
names transformations sparseAssay(),
compactAssay(), disjoinAssay(), and
qreduceAssay(). Because of the complexity of the
transformation, ti usually only makes sense transform
RaggedExperiment objects with a single assay; this is
currently enforced at time of coercion.
Usage
sparseSummarizedExperiment(x, i = 1, withDimnames = TRUE, sparse = FALSE)
compactSummarizedExperiment(x, i = 1L, withDimnames = TRUE, sparse = FALSE)
disjoinSummarizedExperiment(x, simplifyDisjoin, i = 1L, withDimnames = TRUE)
qreduceSummarizedExperiment(
x,
query,
simplifyReduce,
i = 1L,
withDimnames = TRUE
)Arguments
- x
RaggedExperiment- i
integer(1),character(1), orlogical()selecting the assay to be transformed.- withDimnames
logical(1)default TRUE. propagate dimnames to SummarizedExperiment.- sparse
logical(1) whether to return a
sparseMatrixrepresentation- simplifyDisjoin
functionof 1 argument, used to transform assays. Seeassay-functions.- query
GRangesprovding regions over which reduction is to occur.- simplifyReduce
functionof 3 arguments used to transform assays. Seeassay-functions.
Value
All functions return RangedSummarizedExperiment.
sparseSummarizedExperiment has rowRanges()
identical to the row ranges of x, and assay()
data as sparseAssay(). This is very space-inefficient
representation of ragged data. Use 'sparse=TRUE' to obtain
a sparseMatrix assay representation.
compactSummarizedExperiment has rowRanges()
identical to the row ranges of x, and assay()
data as compactAssay(). This is space-inefficient
representation of ragged data when samples are primarily
composed of different ranges. Use 'sparse=TRUE' to obtain
a sparseMatrix assay representation.
disjoinSummarizedExperiment has rowRanges()
identical to the disjoint row ranges of x,
disjoint(rowRanges(x)), and assay() data as
disjoinAssay().
qreduceSummarizedExperiment has rowRanges()
identical to query, and assay() data as
qreduceAssay().
sparseMatrix
Convert a dgCMatrix to a RaggedExperiment given that the rownames
are coercible to GRanges.
In the following example, x is a dgCMatrix from the Matrix package.
`as(x, "RaggedExperiment")`Examples
x <- RaggedExperiment(GRangesList(
GRanges(c("A:1-5", "A:4-6", "A:10-15"), score=1:3),
GRanges(c("A:1-5", "B:1-3"), score=4:5)
))
## sparseSummarizedExperiment
sse <- sparseSummarizedExperiment(x)
assay(sse)
#> [,1] [,2]
#> A:1-5 1 NA
#> A:4-6 2 NA
#> A:10-15 3 NA
#> A:1-5 NA 4
#> B:1-3 NA 5
rowRanges(sse)
#> GRanges object with 5 ranges and 0 metadata columns:
#> seqnames ranges strand
#> <Rle> <IRanges> <Rle>
#> A:1-5 A 1-5 *
#> A:4-6 A 4-6 *
#> A:10-15 A 10-15 *
#> A:1-5 A 1-5 *
#> B:1-3 B 1-3 *
#> -------
#> seqinfo: 2 sequences from an unspecified genome; no seqlengths
## compactSummarizedExperiment
cse <- compactSummarizedExperiment(x)
assay(cse)
#> [,1] [,2]
#> A:1-5 1 4
#> A:4-6 2 NA
#> A:10-15 3 NA
#> B:1-3 NA 5
rowRanges(cse)
#> GRanges object with 4 ranges and 0 metadata columns:
#> seqnames ranges strand
#> <Rle> <IRanges> <Rle>
#> A:1-5 A 1-5 *
#> A:4-6 A 4-6 *
#> A:10-15 A 10-15 *
#> B:1-3 B 1-3 *
#> -------
#> seqinfo: 2 sequences from an unspecified genome; no seqlengths
## disjoinSummarizedExperiment
disjoinAssay(x, lengths)
#> [,1] [,2]
#> A:1-3 1 1
#> A:4-5 2 1
#> A:6 1 NA
#> A:10-15 1 NA
#> B:1-3 NA 1
dse <- disjoinSummarizedExperiment(x, lengths)
assay(dse)
#> [,1] [,2]
#> A:1-3 1 1
#> A:4-5 2 1
#> A:6 1 NA
#> A:10-15 1 NA
#> B:1-3 NA 1
rowRanges(dse)
#> GRanges object with 5 ranges and 0 metadata columns:
#> seqnames ranges strand
#> <Rle> <IRanges> <Rle>
#> A:1-3 A 1-3 *
#> A:4-5 A 4-5 *
#> A:6 A 6 *
#> A:10-15 A 10-15 *
#> B:1-3 B 1-3 *
#> -------
#> seqinfo: 2 sequences from an unspecified genome; no seqlengths
## qreduceSummarizedExperiment
x <- RaggedExperiment(GRangesList(
GRanges(c("A:1-3", "A:4-5", "A:10-15"), score=1:3),
GRanges(c("A:4-5", "B:1-3"), score=4:5)
))
query <- GRanges(c("A:1-2", "A:4-5", "B:1-5"))
weightedmean <- function(scores, ranges, qranges)
{
## weighted average score per query range
## the weight corresponds to the size of the overlap of each
## overlapping subject range with the corresponding query range
isects <- pintersect(ranges, qranges)
sum(scores * width(isects)) / sum(width(isects))
}
qreduceAssay(x, query, weightedmean)
#> [,1] [,2]
#> A:1-2 1 NA
#> A:4-5 2 4
#> B:1-5 NA 5
qse <- qreduceSummarizedExperiment(x, query, weightedmean)
assay(qse)
#> [,1] [,2]
#> A:1-2 1 NA
#> A:4-5 2 4
#> B:1-5 NA 5
rowRanges(qse)
#> GRanges object with 3 ranges and 0 metadata columns:
#> seqnames ranges strand
#> <Rle> <IRanges> <Rle>
#> A:1-2 A 1-2 *
#> A:4-5 A 4-5 *
#> B:1-5 B 1-5 *
#> -------
#> seqinfo: 2 sequences from an unspecified genome; no seqlengths
sm <- Matrix::sparseMatrix(
i = c(2, 3, 4, 3, 4, 3, 4),
j = c(1, 1, 1, 3, 3, 4, 4),
x = c(2L, 4L, 2L, 2L, 2L, 4L, 2L),
dims = c(4, 4),
dimnames = list(
c("chr2:1-10", "chr2:2-10", "chr2:3-10", "chr2:4-10"),
LETTERS[1:4]
)
)
as(sm, "RaggedExperiment")
#> class: RaggedExperiment
#> dim: 7 3
#> assays(1): counts
#> rownames: NULL
#> colnames(3): A C D
#> colData names(0):