Write array blocks
write_block.Rd
Use write_block
to write a block of data to an array-like object.
Note that this function is typically used in the context of block processing
of on-disk objects (e.g. DelayedArray objects), often
in combination with read_block
.
Arguments
- sink
A **writable** array-like object. This is typically a RealizationSink derivative (RealizationSink is a virtual class defined in the DelayedArray package), but not necessarily. See
?RealizationSink
in the DelayedArray package for more information about RealizationSink objects.Although
write_block()
will typically be used on a RealizationSink derivative, it can also be used on an ordinary array or other in-memory array-like object that supports subassignment (`[<-`
), like a SparseArray object from the SparseArray package, or a dgCMatrix object from the Matrix package.- viewport
An ArrayViewport object compatible with
sink
, that is, such thatrefdim(viewport)
is identical todim(sink)
.- block
An array-like object of the same dimensions as
viewport
.
See also
ArrayGrid for ArrayGrid and ArrayViewport objects.
SparseArray objects implemented in the SparseArray package.
read_block
to read a block of data from an array-like object.blockApply
and family, in the DelayedArray package, for convenient block processing of an array-like object.RealizationSink objects implemented in the DelayedArray package for more realistic
write_block
examples.
Examples
## Please note that, although educative, the examples below are somewhat
## artificial and do not illustrate real-world usage of write_block().
## See '?RealizationSink' in the DelayedArray package for more realistic
## read_block/write_block examples.
## ---------------------------------------------------------------------
## BASIC EXAMPLE 1: WRITE A BLOCK TO AN ORDINARY MATRIX (DENSE)
## ---------------------------------------------------------------------
m1 <- matrix(1:30, ncol=5)
m1
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 1 7 13 19 25
#> [2,] 2 8 14 20 26
#> [3,] 3 9 15 21 27
#> [4,] 4 10 16 22 28
#> [5,] 5 11 17 23 29
#> [6,] 6 12 18 24 30
## Define the viewport on 'm1' to write the data to:
block1_dim <- c(4, 3)
viewport1 <- ArrayViewport(dim(m1), IRanges(c(3, 2), width=block1_dim))
viewport1
#> 4 x 3 ArrayViewport object on a 6 x 5 array: [3-6,2-4]
## Data to write:
block1 <- read_block(m1, viewport1) + 1000L
## Write the block:
m1A <- write_block(m1, viewport1, block1)
m1A
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 1 7 13 19 25
#> [2,] 2 8 14 20 26
#> [3,] 3 1009 1015 1021 27
#> [4,] 4 1010 1016 1022 28
#> [5,] 5 1011 1017 1023 29
#> [6,] 6 1012 1018 1024 30
## Sanity checks:
stopifnot(identical(`[<-`(m1, 3:6, 2:4, value=block1), m1A))
m1B <- write_block(m1, viewport1, as(block1, "dgCMatrix"))
stopifnot(identical(m1A, m1B))
## ---------------------------------------------------------------------
## BASIC EXAMPLE 2: WRITE A BLOCK TO A SPARSE MATRIX
## ---------------------------------------------------------------------
m2 <- rsparsematrix(12, 20, density=0.2,
rand.x=function(n) sample(25, n, replace=TRUE))
m2
#> 12 x 20 sparse Matrix of class "dgCMatrix"
#>
#> [1,] . . 24 . . . . . . 24 . 14 . . . 7 11 . . .
#> [2,] . . . 23 . . . . . . . . . . . . 15 . . .
#> [3,] . . 15 . . 2 . . . . 19 . . . . . . . 24 .
#> [4,] . . . 13 . . . . . . . 13 . . . . 15 . . .
#> [5,] . . 14 22 . . . . . . . . . . . . 2 . 19 .
#> [6,] . 7 15 . 19 . 5 . . . . . . . . . . . . .
#> [7,] . 17 . . . 25 . . . . 18 . . . . 15 . . 4 .
#> [8,] . 24 . . . . . . . . . . . . . . . . 5 .
#> [9,] . 18 3 . . 8 . . . . . . . . 24 . . . . .
#> [10,] . 20 9 . . . . 4 . . 16 . . . . . 20 . . .
#> [11,] . . . . . . . 18 . 20 22 . . 3 . 22 . . 20 .
#> [12,] 12 . . 1 . . 21 . . . . . . . . 11 . . . .
## Define the viewport on 'm2' to write the data to:
block2_dim <- c(2, 20)
viewport2 <- ArrayViewport(dim(m2), IRanges(c(1, 1), width=block2_dim))
viewport2
#> 2 x 20 ArrayViewport object on a 12 x 20 array: [1-2, ]
## Data to write:
block2 <- matrix(1001:1040, nrow=2)
## Write the block:
m2A <- write_block(m2, viewport2, block2)
m2A
#> 12 x 20 sparse Matrix of class "dgCMatrix"
#>
#> [1,] 1001 1003 1005 1007 1009 1011 1013 1015 1017 1019 1021 1023 1025 1027
#> [2,] 1002 1004 1006 1008 1010 1012 1014 1016 1018 1020 1022 1024 1026 1028
#> [3,] . . 15 . . 2 . . . . 19 . . .
#> [4,] . . . 13 . . . . . . . 13 . .
#> [5,] . . 14 22 . . . . . . . . . .
#> [6,] . 7 15 . 19 . 5 . . . . . . .
#> [7,] . 17 . . . 25 . . . . 18 . . .
#> [8,] . 24 . . . . . . . . . . . .
#> [9,] . 18 3 . . 8 . . . . . . . .
#> [10,] . 20 9 . . . . 4 . . 16 . . .
#> [11,] . . . . . . . 18 . 20 22 . . 3
#> [12,] 12 . . 1 . . 21 . . . . . . .
#>
#> [1,] 1029 1031 1033 1035 1037 1039
#> [2,] 1030 1032 1034 1036 1038 1040
#> [3,] . . . . 24 .
#> [4,] . . 15 . . .
#> [5,] . . 2 . 19 .
#> [6,] . . . . . .
#> [7,] . 15 . . 4 .
#> [8,] . . . . 5 .
#> [9,] 24 . . . . .
#> [10,] . . 20 . . .
#> [11,] . 22 . . 20 .
#> [12,] . 11 . . . .
## Sanity checks:
stopifnot(identical(`[<-`(m2, 1:2, , value=block2), m2A))
m2B <- write_block(m2, viewport2, as(block2, "dgCMatrix"))
stopifnot(identical(m2A, m2B))
## ---------------------------------------------------------------------
## BASIC EXAMPLE 3: WRITE A BLOCK TO A 3D ARRAY
## ---------------------------------------------------------------------
a3 <- array(1:60, dim=5:3)
## Define the viewport on 'a3' to write the data to:
block3_dim <- c(2, 4, 1)
viewport3 <- ArrayViewport(dim(a3), IRanges(c(1, 1, 3), width=block3_dim))
viewport3
#> 2 x 4 x 1 ArrayViewport object on a 5 x 4 x 3 array: [1-2, ,3]
## Data to write:
block3 <- array(-(1:8), dim=block3_dim)
## Write the block:
a3A <- write_block(a3, viewport3, block3)
a3A
#> , , 1
#>
#> [,1] [,2] [,3] [,4]
#> [1,] 1 6 11 16
#> [2,] 2 7 12 17
#> [3,] 3 8 13 18
#> [4,] 4 9 14 19
#> [5,] 5 10 15 20
#>
#> , , 2
#>
#> [,1] [,2] [,3] [,4]
#> [1,] 21 26 31 36
#> [2,] 22 27 32 37
#> [3,] 23 28 33 38
#> [4,] 24 29 34 39
#> [5,] 25 30 35 40
#>
#> , , 3
#>
#> [,1] [,2] [,3] [,4]
#> [1,] -1 -3 -5 -7
#> [2,] -2 -4 -6 -8
#> [3,] 43 48 53 58
#> [4,] 44 49 54 59
#> [5,] 45 50 55 60
#>
## Sanity checks:
stopifnot(identical(`[<-`(a3, 1:2, , 3, value=block3), a3A))
a3B <- write_block(a3, viewport3, as(block3, "SparseArray"))
stopifnot(identical(a3A, a3B))
## ---------------------------------------------------------------------
## BASIC EXAMPLE 4: WRITE BLOCKS DEFINED BY A GRID
## ---------------------------------------------------------------------
a4 <- array(NA_real_, dim=6:4)
## Define a grid of 2x3x2 blocks on 'a4':
grid4 <- RegularArrayGrid(dim(a4), spacings=c(2,3,2))
grid4
#> 3 x 2 x 2 RegularArrayGrid object on a 6 x 5 x 4 array:
#> , , 1
#>
#> [,1] [,2]
#> [1,] [1-2,1-3,1-2] [1-2,4-5,1-2]
#> [2,] [3-4,1-3,1-2] [3-4,4-5,1-2]
#> [3,] [5-6,1-3,1-2] [5-6,4-5,1-2]
#>
#> , , 2
#>
#> [,1] [,2]
#> [1,] [1-2,1-3,3-4] [1-2,4-5,3-4]
#> [2,] [3-4,1-3,3-4] [3-4,4-5,3-4]
#> [3,] [5-6,1-3,3-4] [5-6,4-5,3-4]
#>
nblock <- length(grid4) # number of blocks
## Walk on the grid and write blocks of random data:
for (bid in seq_len(nblock)) {
viewport <- grid4[[bid]]
block <- array(runif(length(viewport)), dim=dim(viewport))
cat("====== Write block ", bid, "/", nblock, " ======\n", sep="")
a4 <- write_block(a4, viewport, block)
}
#> ====== Write block 1/12 ======
#> ====== Write block 2/12 ======
#> ====== Write block 3/12 ======
#> ====== Write block 4/12 ======
#> ====== Write block 5/12 ======
#> ====== Write block 6/12 ======
#> ====== Write block 7/12 ======
#> ====== Write block 8/12 ======
#> ====== Write block 9/12 ======
#> ====== Write block 10/12 ======
#> ====== Write block 11/12 ======
#> ====== Write block 12/12 ======
a4
#> , , 1
#>
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.3622488 0.6587349 0.2861732 0.5377132 0.08408194
#> [2,] 0.8899170 0.3373870 0.7009997 0.2636477 0.59738269
#> [3,] 0.1440962 0.3997907 0.6791321 0.1465283 0.99485625
#> [4,] 0.5921824 0.6373487 0.9615271 0.9157969 0.99003725
#> [5,] 0.6215599 0.9962952 0.5050395 0.7008730 0.76649118
#> [6,] 0.2452635 0.5527719 0.8481838 0.9370410 0.24177049
#>
#> , , 2
#>
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.8633101 0.29891410 0.7156575 0.22351042 0.6440350
#> [2,] 0.3197115 0.01725763 0.5312104 0.99057988 0.2746681
#> [3,] 0.3395018 0.98237240 0.7998833 0.01387054 0.6368964
#> [4,] 0.7954161 0.27032257 0.6340729 0.80949016 0.9202420
#> [5,] 0.9535708 0.43369423 0.8523407 0.25401791 0.1739656
#> [6,] 0.5334745 0.74230031 0.8391222 0.54708727 0.3132273
#>
#> , , 3
#>
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.6923635 0.5652740 0.4644538 0.008876318 0.14993487
#> [2,] 0.8294130 0.1707580 0.9293323 0.992706620 0.96984320
#> [3,] 0.5831012 0.6734732 0.6573837 0.149600815 0.09887899
#> [4,] 0.6319587 0.7304163 0.5192295 0.149469787 0.97117306
#> [5,] 0.6778980 0.2099045 0.1583046 0.606963684 0.96536046
#> [6,] 0.9982615 0.1727390 0.3540038 0.307511342 0.31445338
#>
#> , , 4
#>
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.9924558 0.8719124 0.7872488 0.81415289 0.6236400
#> [2,] 0.9999725 0.6038597 0.6038978 0.04114619 0.7031565
#> [3,] 0.2104369 0.1479698 0.4605602 0.55070559 0.9637677
#> [4,] 0.7931375 0.5087944 0.6155133 0.83027343 0.1572618
#> [5,] 0.9736047 0.8224273 0.4644698 0.41317681 0.1073284
#> [6,] 0.8178529 0.4098893 0.6270261 0.77479378 0.6880921
#>
## ---------------------------------------------------------------------
## BASIC EXAMPLE 5: READ, PROCESS, AND WRITE BLOCKS DEFINED BY TWO GRIDS
## ---------------------------------------------------------------------
## Say we have a 3D array and want to collapse its 3rd dimension by
## summing the array elements that are stacked vertically, that is, we
## want to compute the matrix 'm' obtained by doing 'sum(a[i, j, ])' for
## all valid i and j. There are several ways to do this.
## 1. Here is a solution based on apply():
collapse_3rd_dim <- function(a) apply(a, MARGIN=1:2, sum)
## 2. Here is a slightly more efficient solution:
collapse_3rd_dim <- function(a) {
m <- matrix(0, nrow=nrow(a), ncol=ncol(a))
for (z in seq_len(dim(a)[[3]]))
m <- m + a[ , , z]
m
}
## 3. And here is a block-processing solution that involves two grids,
## one for the sink, and one for the input:
a5 <- array(runif(8000), dim=c(25, 40, 8)) # input
m <- array(NA_real_, dim=dim(a5)[1:2]) # sink
## Since we're going to walk on the two grids simultaneously, read a
## block from 'a5' and write it to 'm', we need to make sure that we
## define grids that are "aligned". More precisely, the two grids must
## have the same number of viewports, and the viewports in one must
## correspond to the viewports in the other one:
m_grid <- RegularArrayGrid(dim(m), spacings=c(10, 10))
a5_grid <- RegularArrayGrid(dim(a5), spacings=c(10, 10, dim(a5)[[3]]))
## Let's check that our two grids are actually "aligned":
stopifnot(identical(length(m_grid), length(a5_grid)))
stopifnot(identical(dims(m_grid), dims(a5_grid)[ , 1:2, drop=FALSE]))
## Walk on the two grids simultaneously, and read/collapse/write blocks:
for (bid in seq_along(m_grid)) {
## Read block from 'a5'.
a5_viewport <- a5_grid[[bid]]
block <- read_block(a5, a5_viewport)
## Collapse it.
block <- collapse_3rd_dim(block)
## Write the collapsed block to 'm'.
m_viewport <- m_grid[[bid]]
m <- write_block(m, m_viewport, block)
}
## Sanity checks:
stopifnot(identical(dim(a5)[1:2], dim(m)))
stopifnot(identical(sum(a5), sum(m)))
stopifnot(identical(collapse_3rd_dim(a5), m))
## See '?RealizationSink' in the DelayedArray package for a more
## realistic "array collapse" example where the blocks are written
## to a RealizationSink object.