DelayedArray objects
DelayedArray-class.RdWrapping an array-like object (typically an on-disk object) in a DelayedArray object allows one to perform common array operations on it without loading the object in memory. In order to reduce memory usage and optimize performance, operations on the object are either delayed or executed using a block processing mechanism.
In-memory versus on-disk realization
To realize a DelayedArray object (i.e. to trigger execution of the
delayed operations carried by the object and return the result as an
ordinary array), call as.array on it. However this realizes the
full object at once in memory which could require too much memory
if the object is big. A big DelayedArray object is preferrably realized
on disk e.g. by calling writeHDF5Array on
it (this function is defined in the HDF5Array package) or coercing it
to an HDF5Array object with as(x, "HDF5Array").
Other on-disk backends can be supported. This uses a block processing
strategy so that the full object is not realized at once in memory. Instead
the object is processed block by block i.e. the blocks are realized in
memory and written to disk one at a time.
See ?writeHDF5Array in the HDF5Array package
for more information about this.
Accessors
DelayedArray objects support the same set of getters as ordinary arrays
i.e. dim(), length(), and dimnames().
In addition, they support type(), nseed(),
seed(), and path().
type() is the DelayedArray equivalent of typeof() (or
storage.mode()) for ordinary arrays and vectors. Note that, for
convenience and consistency, type() also supports ordinary arrays
and vectors. It should also support any array-like object, that is, any
object x for which dim(x) is not NULL.
dimnames(), seed(), and path() also
work as setters.
Subsetting
A DelayedArray object can be subsetted with [ like an ordinary array,
but with the following differences:
N-dimensional single bracket subsetting (i.e. subsetting of the form
x[i_1, i_2, ..., i_n]with one (possibly missing) subscript per dimension) returns a DelayedArray object where the subsetting is actually delayed. So it's a very light operation. One notable exception is whendrop=TRUEand the result has only one dimension, in which case it is realized as an ordinary vector (atomic or list). Note that NAs in the subscripts are not supported.1D-style single bracket subsetting (i.e. subsetting of the form
x[i]) only works if the subscriptiis a numeric or logical vector, or a logical array-like object with the same dimensions asx, or a numeric matrix with one column per dimension inx. Wheniis a numeric vector, all the indices in it must be >= 1 and <=length(x). NAs in the subscripts are not supported. This is NOT a delayed operation (block processing is triggered) i.e. the result is realized as an ordinary vector (atomic or list). One exception is whenxhas only one dimension anddropis set toFALSE, in which case the subsetting is delayed.
Subsetting with [[ is supported but only the 1D-style form of it
at the moment, that is, subsetting of the form x[[i]] where i
is a single numeric value >= 1 and <= length(x). It is
equivalent to x[i][[1]].
Subassignment to a DelayedArray object with [<- is also supported
like with an ordinary array, but with the following restrictions:
N-dimensional subassignment (i.e. subassignment of the form
x[i_1, i_2, ..., i_n] <- valuewith one (possibly missing) subscript per dimension) only accepts a replacement value (a.k.a. right value) that is an array-like object (e.g. ordinary array, dgCMatrix object, DelayedArray object, etc...) or an ordinary vector (atomic or list) of length 1.1D-style subassignment (a.k.a. 1D-style subassignment, that is, subassignment of the form
x[i] <- value) only works if the subscriptiis a logical DelayedArray object of the same dimensions asxand if the replacement value is an ordinary vector (atomic or list) of length 1.Filling with a vector, that is, subassignment of the form
x[] <- vwherevis an ordinary vector (atomic or list), is only supported if the length of the vector is a divisor ofnrow(x).
These 3 forms of subassignment are implemented as delayed operations so are very light.
Single value replacement (x[[...]] <- value) is not supported yet.
See also
showtreefor DelayedArray accessorsnseed,seed, andpath.realizefor realizing a DelayedArray object in memory or on disk.blockApplyand family for convenient block processing of an array-like object.DelayedArray-utils for common operations on DelayedArray objects.
DelayedArray-stats for statistical functions on DelayedArray objects.
matrixStats-methods for DelayedMatrix row/col summarization.
DelayedMatrix-rowsum for
rowsum()andcolsum()methods for DelayedMatrix objects.DelayedMatrix-mult for DelayedMatrix multiplication and cross-product.
ConstantArray objects for mimicking an array containing a constant value, without actually creating said array in memory.
RleArray objects for representing in-memory Run Length Encoded array-like datasets.
HDF5Array objects in the HDF5Array package.
DataFrame objects in the S4Vectors package.
array objects in base R.
Examples
## ---------------------------------------------------------------------
## A. WRAP AN ORDINARY ARRAY IN A DelayedArray OBJECT
## ---------------------------------------------------------------------
a <- array(runif(1500000), dim=c(10000, 30, 5))
A <- DelayedArray(a)
A
#> <10000 x 30 x 5> DelayedArray object of type "double":
#> ,,1
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] 0.17295454 0.50228111 0.08797905 . 0.9959464 0.8348479
#> [2,] 0.06493924 0.97783460 0.96461208 . 0.9603776 0.5256657
#> ... . . . . . .
#> [9999,] 0.69219573 0.07848326 0.07278009 . 0.3675142 0.6797809
#> [10000,] 0.29900256 0.22960714 0.05297260 . 0.0530765 0.6870121
#>
#> ...
#>
#> ,,5
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] 0.85008834 0.29234868 0.47663892 . 0.2109716 0.4997599
#> [2,] 0.03570222 0.15024342 0.20017364 . 0.8747636 0.7099284
#> ... . . . . . .
#> [9999,] 0.3655142 0.7945577 0.1624435 . 0.8417917 0.7864408
#> [10000,] 0.1596408 0.4313499 0.4138495 . 0.6194634 0.9775811
#>
## The seed of a DelayedArray object is **always** treated as a
## "read-only" object so will never be modified by the operations
## we perform on A:
stopifnot(identical(a, seed(A)))
type(A)
#> [1] "double"
## N-dimensional single bracket subsetting:
m <- a[11:20 , 5, -3] # an ordinary matrix
M <- A[11:20 , 5, -3] # a DelayedMatrix object
stopifnot(identical(m, as.array(M)))
## 1D-style single bracket subsetting:
A[11:20]
#> [1] 0.07758068 0.59965366 0.87901081 0.73203716 0.08005309 0.83274981
#> [7] 0.86470570 0.71296497 0.75740711 0.12243392
A[A <= 1e-5]
#> [1] 4.966976e-06 9.604031e-06 7.921131e-06 1.286156e-06 3.997236e-06
#> [6] 9.523937e-06 6.426126e-08 4.861970e-06 7.504830e-06 2.031447e-06
#> [11] 7.894589e-06 9.104144e-06 2.589077e-06 5.005859e-08 7.156515e-06
#> [16] 9.323005e-06 6.955815e-06 4.800968e-06 2.797926e-06 9.188661e-06
stopifnot(identical(a[a <= 1e-5], A[A <= 1e-5]))
## Subassignment:
A[A < 0.2] <- NA
a[a < 0.2] <- NA
stopifnot(identical(a, as.array(A)))
A[2:5, 1:2, ] <- array(1:40, c(4, 2, 5))
a[2:5, 1:2, ] <- array(1:40, c(4, 2, 5))
stopifnot(identical(a, as.array(A)))
## Other operations:
crazy <- function(x) (5 * x[ , , 1] ^ 3 + 1L) * log(x[, , 2])
b <- crazy(a)
head(b)
#> [,1] [,2] [,3] [,4] [,5] [,6]
#> [1,] NA -0.23948556 NA -3.20707020 NA -1.1470040
#> [2,] 13.18335 1605.65829777 -6.0060106 -0.11856540 NA -1.5028617
#> [3,] 94.40599 2852.82097331 NA -0.01366113 NA -0.2739800
#> [4,] 326.11376 4647.01414509 NA NA -3.0554730 NA
#> [5,] 797.65503 7100.59971766 -0.2557321 -5.39698626 NA -2.2827165
#> [6,] NA -0.01075831 -1.6697196 NA -0.7936828 -0.2565307
#> [,7] [,8] [,9] [,10] [,11] [,12]
#> [1,] -3.6022805 -1.9897722 -0.5379044 -0.1123039 NA -1.409749
#> [2,] NA -1.3275820 NA NA -1.5070459 NA
#> [3,] -0.8327830 -0.6767560 NA NA -2.2313879 NA
#> [4,] -0.2822528 -0.5485830 NA -0.1162019 -0.3672875 -1.542087
#> [5,] -0.3206836 -0.6488671 NA -0.9619093 -2.4693303 NA
#> [6,] NA NA -7.8777190 NA -0.4903264 NA
#> [,13] [,14] [,15] [,16] [,17] [,18]
#> [1,] NA NA -0.09224164 NA -1.7857609 -1.3962225
#> [2,] -0.08518998 -0.7431843 NA NA -1.7618634 -0.8066919
#> [3,] NA -0.7269428 NA -0.02858646 -0.6036519 -1.2630364
#> [4,] -6.59402347 NA NA -1.22819253 NA NA
#> [5,] NA NA -2.99159355 -1.12541116 -0.1268075 -0.2867510
#> [6,] -0.33978357 -0.3737906 NA NA -5.3907519 -0.2840287
#> [,19] [,20] [,21] [,22] [,23] [,24]
#> [1,] NA -2.85748962 -0.6360317 NA NA -0.3039220
#> [2,] NA -0.45447829 -1.5708661 NA -2.6177232 NA
#> [3,] -1.0949013 -2.04682980 -0.6451412 NA -1.1834395 NA
#> [4,] NA -1.96302734 -0.7520166 -1.351614 -0.9731677 -0.3866473
#> [5,] -0.3998895 NA NA -1.588815 -4.6645790 -2.4955097
#> [6,] -1.1495873 -0.05182364 -0.3109204 -2.188424 -0.7635950 -1.2213688
#> [,25] [,26] [,27] [,28] [,29] [,30]
#> [1,] -7.329597 -2.9902101 NA -0.4733355 -0.6566356 -0.714658
#> [2,] -7.397454 -2.4980368 NA NA -1.0505077 -1.281317
#> [3,] NA NA NA NA -0.2307934 NA
#> [4,] NA -0.6705747 NA NA NA -1.466592
#> [5,] NA NA -0.2745992 -7.2925156 NA -1.970912
#> [6,] NA -0.5712061 -1.6221770 NA NA -1.281945
B <- crazy(A) # very fast! (all operations are delayed)
B
#> <10000 x 30> DelayedMatrix object of type "double":
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] NA -0.2394856 NA . -0.6566356 -0.7146580
#> [2,] 13.1833475 1605.6582978 -6.0060106 . -1.0505077 -1.2813167
#> [3,] 94.4059888 2852.8209733 NA . -0.2307934 NA
#> [4,] 326.1137571 4647.0141451 NA . NA -1.4665921
#> [5,] 797.6550346 7100.5997177 -0.2557321 . NA -1.9709123
#> ... . . . . . .
#> [9996,] -0.34780449 -1.70031354 NA . NA NA
#> [9997,] -0.56893931 NA NA . NA NA
#> [9998,] -2.34750677 -1.96438054 NA . NA -0.534062
#> [9999,] -0.02652313 NA NA . NA -3.637590
#> [10000,] -0.47849404 -0.79620916 NA . NA -2.816963
cs <- colSums(b)
CS <- colSums(B)
stopifnot(identical(cs, CS))
## ---------------------------------------------------------------------
## B. WRAP A DataFrame OBJECT IN A DelayedArray OBJECT
## ---------------------------------------------------------------------
## Generate random coverage and score along an imaginary chromosome:
cov <- Rle(sample(20, 5000, replace=TRUE), sample(6, 5000, replace=TRUE))
score <- Rle(sample(100, nrun(cov), replace=TRUE), runLength(cov))
DF <- DataFrame(cov, score)
A2 <- DelayedArray(DF)
A2
#> <17684 x 2> DelayedMatrix object of type "integer":
#> cov score
#> [1,] 2 20
#> [2,] 2 20
#> [3,] 20 52
#> [4,] 20 52
#> [5,] 20 52
#> ... . .
#> [17680,] 10 16
#> [17681,] 10 16
#> [17682,] 10 16
#> [17683,] 10 16
#> [17684,] 10 16
seed(A2) # 'DF'
#> DataFrame with 17684 rows and 2 columns
#> cov score
#> <Rle> <Rle>
#> 1 2 20
#> 2 2 20
#> 3 20 52
#> 4 20 52
#> 5 20 52
#> ... ... ...
#> 17680 10 16
#> 17681 10 16
#> 17682 10 16
#> 17683 10 16
#> 17684 10 16
## Coercion of a DelayedMatrix object to DataFrame produces a DataFrame
## object with Rle columns:
as(A2, "DataFrame")
#> DataFrame with 17684 rows and 2 columns
#> cov score
#> <Rle> <Rle>
#> 1 2 20
#> 2 2 20
#> 3 20 52
#> 4 20 52
#> 5 20 52
#> ... ... ...
#> 17680 10 16
#> 17681 10 16
#> 17682 10 16
#> 17683 10 16
#> 17684 10 16
stopifnot(identical(DF, as(A2, "DataFrame")))
t(A2) # transposition is delayed so is very fast and memory-efficient
#> <2 x 17684> DelayedMatrix object of type "integer":
#> [,1] [,2] [,3] [,4] ... [,17681] [,17682] [,17683]
#> cov 2 2 20 20 . 10 10 10
#> score 20 20 52 52 . 16 16 16
#> [,17684]
#> cov 10
#> score 16
colSums(A2)
#> cov score
#> 185010 879940
## ---------------------------------------------------------------------
## C. AN HDF5Array OBJECT IS A (PARTICULAR KIND OF) DelayedArray OBJECT
## ---------------------------------------------------------------------
library(HDF5Array)
A3 <- as(a, "HDF5Array") # write 'a' to an HDF5 file
A3
#> <10000 x 30 x 5> HDF5Array object of type "double":
#> ,,1
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] NA 0.5022811 NA . 0.9959464 0.8348479
#> [2,] 1.0000000 5.0000000 0.9646121 . 0.9603776 0.5256657
#> ... . . . . . .
#> [9999,] 0.6921957 NA NA . 0.3675142 0.6797809
#> [10000,] 0.2990026 0.2296071 NA . NA 0.6870121
#>
#> ...
#>
#> ,,5
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] 0.8500883 0.2923487 0.4766389 . 0.2109716 0.4997599
#> [2,] 33.0000000 37.0000000 0.2001736 . 0.8747636 0.7099284
#> ... . . . . . .
#> [9999,] 0.3655142 0.7945577 NA . 0.8417917 0.7864408
#> [10000,] NA 0.4313499 0.4138495 . 0.6194634 0.9775811
#>
is(A3, "DelayedArray") # TRUE
#> [1] TRUE
seed(A3) # an HDF5ArraySeed object
#> An object of class "HDF5ArraySeed"
#> Slot "filepath":
#> [1] "/tmp/RtmptRnwDm/HDF5Array_dump/auto1ce158ac37d3.h5"
#>
#> Slot "name":
#> [1] "/HDF5ArrayAUTO00002"
#>
#> Slot "as_sparse":
#> [1] FALSE
#>
#> Slot "type":
#> [1] NA
#>
#> Slot "dim":
#> [1] 10000 30 5
#>
#> Slot "chunkdim":
#> [1] 8735 26 4
#>
#> Slot "first_val":
#> [1] NA
#>
B3 <- crazy(A3) # very fast! (all operations are delayed)
B3 # not an HDF5Array object anymore because
#> <10000 x 30> DelayedMatrix object of type "double":
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] NA -0.2394856 NA . -0.6566356 -0.7146580
#> [2,] 13.1833475 1605.6582978 -6.0060106 . -1.0505077 -1.2813167
#> [3,] 94.4059888 2852.8209733 NA . -0.2307934 NA
#> [4,] 326.1137571 4647.0141451 NA . NA -1.4665921
#> [5,] 797.6550346 7100.5997177 -0.2557321 . NA -1.9709123
#> ... . . . . . .
#> [9996,] -0.34780449 -1.70031354 NA . NA NA
#> [9997,] -0.56893931 NA NA . NA NA
#> [9998,] -2.34750677 -1.96438054 NA . NA -0.534062
#> [9999,] -0.02652313 NA NA . NA -3.637590
#> [10000,] -0.47849404 -0.79620916 NA . NA -2.816963
# now it carries delayed operations
CS3 <- colSums(B3)
stopifnot(identical(cs, CS3))
## ---------------------------------------------------------------------
## D. PERFORM THE DELAYED OPERATIONS
## ---------------------------------------------------------------------
as(B3, "HDF5Array") # "realize" 'B3' on disk
#> <10000 x 30> HDF5Matrix object of type "double":
#> [,1] [,2] [,3] ... [,29] [,30]
#> [1,] NA -0.2394856 NA . -0.6566356 -0.7146580
#> [2,] 13.1833475 1605.6582978 -6.0060106 . -1.0505077 -1.2813167
#> [3,] 94.4059888 2852.8209733 NA . -0.2307934 NA
#> [4,] 326.1137571 4647.0141451 NA . NA -1.4665921
#> [5,] 797.6550346 7100.5997177 -0.2557321 . NA -1.9709123
#> ... . . . . . .
#> [9996,] -0.34780449 -1.70031354 NA . NA NA
#> [9997,] -0.56893931 NA NA . NA NA
#> [9998,] -2.34750677 -1.96438054 NA . NA -0.534062
#> [9999,] -0.02652313 NA NA . NA -3.637590
#> [10000,] -0.47849404 -0.79620916 NA . NA -2.816963
## If this is just an intermediate result, you can either keep going
## with B3 or replace it with its "realized" version:
B3 <- as(B3, "HDF5Array") # no more delayed operations on new 'B3'
seed(B3)
#> An object of class "HDF5ArraySeed"
#> Slot "filepath":
#> [1] "/tmp/RtmptRnwDm/HDF5Array_dump/auto1ce11f034507.h5"
#>
#> Slot "name":
#> [1] "/HDF5ArrayAUTO00004"
#>
#> Slot "as_sparse":
#> [1] FALSE
#>
#> Slot "type":
#> [1] NA
#>
#> Slot "dim":
#> [1] 10000 30
#>
#> Slot "chunkdim":
#> [1] 10000 30
#>
#> Slot "first_val":
#> [1] NA
#>
path(B3)
#> [1] "/tmp/RtmptRnwDm/HDF5Array_dump/auto1ce11f034507.h5"
## For convenience, realize() can be used instead of explicit coercion.
## The current "automatic realization backend" controls where
## realization happens e.g. in memory if set to NULL or in an HDF5
## file if set to "HDF5Array":
D <- cbind(B3, exp(B3))
D
#> <10000 x 60> DelayedMatrix object of type "double":
#> [,1] [,2] [,3] ... [,59] [,60]
#> [1,] NA -0.2394856 NA . 0.5185931 0.4893594
#> [2,] 13.1833475 1605.6582978 -6.0060106 . 0.3497601 0.2776715
#> [3,] 94.4059888 2852.8209733 NA . 0.7939034 NA
#> [4,] 326.1137571 4647.0141451 NA . NA 0.2307104
#> [5,] 797.6550346 7100.5997177 -0.2557321 . NA 0.1393297
#> ... . . . . . .
#> [9996,] -0.34780449 -1.70031354 NA . NA NA
#> [9997,] -0.56893931 NA NA . NA NA
#> [9998,] -2.34750677 -1.96438054 NA . NA 0.58621889
#> [9999,] -0.02652313 NA NA . NA 0.02631569
#> [10000,] -0.47849404 -0.79620916 NA . NA 0.05978724
setAutoRealizationBackend("HDF5Array")
D <- realize(D)
D
#> <10000 x 60> HDF5Matrix object of type "double":
#> [,1] [,2] [,3] ... [,59] [,60]
#> [1,] NA -0.2394856 NA . 0.5185931 0.4893594
#> [2,] 13.1833475 1605.6582978 -6.0060106 . 0.3497601 0.2776715
#> [3,] 94.4059888 2852.8209733 NA . 0.7939034 NA
#> [4,] 326.1137571 4647.0141451 NA . NA 0.2307104
#> [5,] 797.6550346 7100.5997177 -0.2557321 . NA 0.1393297
#> ... . . . . . .
#> [9996,] -0.34780449 -1.70031354 NA . NA NA
#> [9997,] -0.56893931 NA NA . NA NA
#> [9998,] -2.34750677 -1.96438054 NA . NA 0.58621889
#> [9999,] -0.02652313 NA NA . NA 0.02631569
#> [10000,] -0.47849404 -0.79620916 NA . NA 0.05978724
## See '?setAutoRealizationBackend' for more information about
## "realization backends".
setAutoRealizationBackend() # restore default (NULL)
## ---------------------------------------------------------------------
## E. MODIFY THE PATH OF A DelayedArray OBJECT
## ---------------------------------------------------------------------
## This can be useful if the file containing the array data is on a
## shared partition but the exact path to the partition depends on the
## machine from which the data is being accessed.
## For example:
if (FALSE) { # \dontrun{
library(HDF5Array)
A <- HDF5Array("/path/to/lab_data/my_precious_data.h5")
path(A)
## Operate on A...
## Now A carries delayed operations.
## Make sure path(A) still works:
path(A)
## Save A:
save(A, file="A.rda")
## A.rda should be small (it doesn't contain the array data).
## Send it to a co-worker that has access to my_precious_data.h5.
## Co-worker loads it:
load("A.rda")
path(A)
## A is broken because path(A) is incorrect for co-worker:
A # error!
## Co-worker fixes the path (in this case this is better done using the
## dirname() setter rather than the path() setter):
dirname(A) <- "E:/other/path/to/lab_data"
## A "works" again:
A
} # }
## ---------------------------------------------------------------------
## F. WRAP A SPARSE MATRIX IN A DelayedArray OBJECT
## ---------------------------------------------------------------------
if (FALSE) { # \dontrun{
M <- 75000L
N <- 1800L
p <- sparseMatrix(sample(M, 9000000, replace=TRUE),
sample(N, 9000000, replace=TRUE),
x=runif(9000000), dims=c(M, N))
P <- DelayedArray(p)
P
p2 <- as(P, "sparseMatrix")
stopifnot(identical(p, p2))
## The following is based on the following post by Murat Tasan on the
## R-help mailing list:
## https://stat.ethz.ch/pipermail/r-help/2017-May/446702.html
## As pointed out by Murat, the straight-forward row normalization
## directly on sparse matrix 'p' would consume too much memory:
row_normalized_p <- p / rowSums(p^2) # consumes too much memory
## because the rowSums() result is being recycled (appropriately) into a
## *dense* matrix with dimensions equal to dim(p).
## Murat came up with the following solution that is very fast and
## memory-efficient:
row_normalized_p1 <- Diagonal(x=1/sqrt(Matrix::rowSums(p^2)))
## With a DelayedArray object, the straight-forward approach uses a
## block processing strategy behind the scene so it doesn't consume
## too much memory.
## First, let's see block processing in action:
DelayedArray:::set_verbose_block_processing(TRUE)
## and check the automatic block size:
getAutoBlockSize()
row_normalized_P <- P / sqrt(DelayedArray::rowSums(P^2))
## Increasing the block size increases the speed but also memory usage:
setAutoBlockSize(2e8)
row_normalized_P2 <- P / sqrt(DelayedArray::rowSums(P^2))
stopifnot(all.equal(row_normalized_P, row_normalized_P2))
## Back to sparse representation:
DelayedArray:::set_verbose_block_processing(FALSE)
row_normalized_p2 <- as(row_normalized_P, "sparseMatrix")
stopifnot(all.equal(row_normalized_p1, row_normalized_p2))
setAutoBlockSize() # reset automatic block size to factory settings
} # }