mat1 <- matrix(1:10, ncol = 2)
mat2 <- matrix(2:11, ncol = 2)
list_matrices <- list(mat1, mat2)
list_matrices
#> [[1]]
#> [,1] [,2]
#> [1,] 1 6
#> [2,] 2 7
#> [3,] 3 8
#> [4,] 4 9
#> [5,] 5 10
#>
#> [[2]]
#> [,1] [,2]
#> [1,] 2 7
#> [2,] 3 8
#> [3,] 4 9
#> [4,] 5 10
#> [5,] 6 11
out <- do.call(`+`, list_matrices) / length(list_matrices)
out
#> [,1] [,2]
#> [1,] 1.5 6.5
#> [2,] 2.5 7.5
#> [3,] 3.5 8.5
#> [4,] 4.5 9.5
#> [5,] 5.5 10.5
# Can use `Reduce` instead of `do.call` ?
Reduce(`+`, list_matrices) / length(list_matrices)
#> [,1] [,2]
#> [1,] 1.5 6.5
#> [2,] 2.5 7.5
#> [3,] 3.5 8.5
#> [4,] 4.5 9.5
#> [5,] 5.5 10.5
# you can also instead construct the matrices as arras and then do rowMeans
new_arrays <- simplify2array(list_matrices)
new_arrays
#> , , 1
#>
#> [,1] [,2]
#> [1,] 1 6
#> [2,] 2 7
#> [3,] 3 8
#> [4,] 4 9
#> [5,] 5 10
#>
#> , , 2
#>
#> [,1] [,2]
#> [1,] 2 7
#> [2,] 3 8
#> [3,] 4 9
#> [4,] 5 10
#> [5,] 6 11
rowMeans(x = new_arrays, dims = 2)
#> [,1] [,2]
#> [1,] 1.5 6.5
#> [2,] 2.5 7.5
#> [3,] 3.5 8.5
#> [4,] 4.5 9.5
#> [5,] 5.5 10.5
# for completeness, here's a benchmark
mean_do_call <- function(list_matrices) {
do.call(`+`, list_matrices) / length(list_matrices)
}
mean_reduce <- function(list_matrices) {
Reduce(`+`, list_matrices) / length(list_matrices)
}
mean_rowMeans <- function(mat_array) {
rowMeans(mat_array, dims = 2)
}
bm <- bench::mark(
do_call = mean_do_call(list_matrices),
reduce = mean_reduce(list_matrices),
rowmeans = mean_rowMeans(new_arrays)
)
bm
#> # A tibble: 3 × 6
#> expression min median `itr/sec` mem_alloc `gc/sec`
#> <bch:expr> <bch:tm> <bch:tm> <dbl> <bch:byt> <dbl>
#> 1 do_call 983.94ns 1.15µs 727550. 0B 72.8
#> 2 reduce 1.72µs 1.97µs 471671. 0B 47.2
#> 3 rowmeans 1.84µs 2.09µs 440274. 0B 44.0
summary(bm, relative = TRUE)
#> # A tibble: 3 × 6
#> expression min median `itr/sec` mem_alloc `gc/sec`
#> <bch:expr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 do_call 1 1 1.65 NaN 1.65
#> 2 reduce 1.75 1.71 1.07 NaN 1.07
#> 3 rowmeans 1.88 1.82 1 NaN 1
plot(bm)
#> Loading required namespace: tidyr
# what about for larger matrices?
bigger_matrix <- function(square_size) {
n_values <- square_size * square_size
vals <- runif(n_values)
matrix(data = vals, nrow = square_size, ncol = square_size)
}
b_mat1 <- bigger_matrix(1000)
b_mat2 <- bigger_matrix(1000)
dim(b_mat1)
#> [1] 1000 1000
b_list_matrices <- list(b_mat1, b_mat2)
b_new_arrays <- simplify2array(b_list_matrices)
dim(b_new_arrays)
#> [1] 1000 1000 2
big_bm <- bench::mark(
do_call = mean_do_call(b_list_matrices),
reduce = mean_reduce(b_list_matrices),
rowmeans = mean_rowMeans(b_new_arrays)
)
big_bm
#> # A tibble: 3 × 6
#> expression min median `itr/sec` mem_alloc `gc/sec`
#> <bch:expr> <bch:tm> <bch:tm> <dbl> <bch:byt> <dbl>
#> 1 do_call 1.55ms 2.01ms 386. 7.63MB 102.
#> 2 reduce 1.52ms 1.94ms 485. 7.63MB 70.1
#> 3 rowmeans 1.68ms 2.5ms 301. 7.63MB 45.0
summary(big_bm, relative = TRUE)
#> # A tibble: 3 × 6
#> expression min median `itr/sec` mem_alloc `gc/sec`
#> <bch:expr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 do_call 1.02 1.04 1.28 1 2.27
#> 2 reduce 1 1 1.61 1 1.56
#> 3 rowmeans 1.10 1.29 1 1 1
plot(big_bm)
# let's explore this across a range of values using bench press, just because
bp <- bench::press(
square_size = c(10, 100, seq(from = 1000, to = 5000, by = 1000)),
reps = 1,
{
b_mat1 <- bigger_matrix(square_size)
b_mat2 <- bigger_matrix(square_size)
b_list_matrices <- list(b_mat1, b_mat2)
b_new_arrays <- simplify2array(b_list_matrices)
bench::mark(
do_call = mean_do_call(b_list_matrices),
reduce = mean_reduce(b_list_matrices),
rowmeans = mean_rowMeans(b_new_arrays)
)
}
)
#> Running with:
#> square_size reps
#> 1 10 1
#> 2 100 1
#> 3 1000 1
#> 4 2000 1
#> 5 3000 1
#> 6 4000 1
#> 7 5000 1
bp
#> # A tibble: 21 × 8
#> expression square_size reps min median `itr/sec` mem_alloc `gc/sec`
#> <bch:expr> <dbl> <dbl> <bch:tm> <bch:tm> <dbl> <bch:byt> <dbl>
#> 1 do_call 10 1 1.11µs 1.23µs 735988. 848B 0
#> 2 reduce 10 1 1.89µs 2.09µs 456900. 848B 45.7
#> 3 rowmeans 10 1 1.97µs 2.3µs 316825. 848B 0
#> 4 do_call 100 1 14.14µs 19.76µs 42885. 78.17KB 51.5
#> 5 reduce 100 1 14.8µs 19.76µs 50505. 78.17KB 60.7
#> 6 rowmeans 100 1 12.55µs 16.89µs 58505. 78.17KB 76.2
#> 7 do_call 1000 1 1.52ms 1.93ms 497. 7.63MB 45.1
#> 8 reduce 1000 1 1.51ms 2.04ms 481. 7.63MB 44.6
#> 9 rowmeans 1000 1 1.22ms 2.1ms 477. 7.63MB 44.9
#> 10 do_call 2000 1 7.62ms 7.77ms 124. 30.52MB 28.9
#> # ℹ 11 more rows
plot(bp)
Created on 2025-03-14 with reprex v2.1.1
Session info
sessioninfo::session_info()
#> ─ Session info ───────────────────────────────────────────────────────────────
#> setting value
#> version R version 4.4.2 (2024-10-31)
#> os macOS Sequoia 15.1
#> system aarch64, darwin20
#> ui X11
#> language (EN)
#> collate en_US.UTF-8
#> ctype en_US.UTF-8
#> tz Australia/Hobart
#> date 2025-03-14
#> pandoc 3.2.1 @ /opt/homebrew/bin/ (via rmarkdown)
#>
#> ─ Packages ───────────────────────────────────────────────────────────────────
#> package * version date (UTC) lib source
#> beeswarm 0.4.0 2021-06-01 [1] CRAN (R 4.4.0)
#> bench 1.1.3 2023-05-04 [1] CRAN (R 4.4.0)
#> cli 3.6.4 2025-02-13 [1] CRAN (R 4.4.1)
#> colorspace 2.1-1 2024-07-26 [1] CRAN (R 4.4.0)
#> curl 6.2.0 2025-01-23 [1] CRAN (R 4.4.1)
#> digest 0.6.37 2024-08-19 [1] CRAN (R 4.4.1)
#> dplyr 1.1.4 2023-11-17 [1] CRAN (R 4.4.0)
#> evaluate 1.0.1 2024-10-10 [1] CRAN (R 4.4.1)
#> farver 2.1.2 2024-05-13 [1] CRAN (R 4.4.0)
#> fastmap 1.2.0 2024-05-15 [1] CRAN (R 4.4.0)
#> fs 1.6.5 2024-10-30 [1] CRAN (R 4.4.1)
#> generics 0.1.3 2022-07-05 [1] CRAN (R 4.4.0)
#> ggbeeswarm 0.7.2 2023-04-29 [1] CRAN (R 4.4.0)
#> ggplot2 3.5.1 2024-04-23 [1] CRAN (R 4.4.0)
#> glue 1.8.0 2024-09-30 [1] CRAN (R 4.4.1)
#> gtable 0.3.6 2024-10-25 [1] CRAN (R 4.4.1)
#> htmltools 0.5.8.1 2024-04-04 [1] CRAN (R 4.4.0)
#> knitr 1.49 2024-11-08 [1] CRAN (R 4.4.1)
#> labeling 0.4.3 2023-08-29 [1] CRAN (R 4.4.0)
#> lifecycle 1.0.4 2023-11-07 [1] CRAN (R 4.4.0)
#> magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.4.0)
#> munsell 0.5.1 2024-04-01 [1] CRAN (R 4.4.0)
#> pillar 1.10.1 2025-01-07 [1] CRAN (R 4.4.1)
#> pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.4.0)
#> profmem 0.6.0 2020-12-13 [1] CRAN (R 4.4.0)
#> purrr 1.0.4 2025-02-05 [1] CRAN (R 4.4.1)
#> R6 2.6.1 2025-02-15 [1] CRAN (R 4.4.1)
#> reprex 2.1.1 2024-07-06 [1] CRAN (R 4.4.0)
#> rlang 1.1.5 2025-01-17 [1] CRAN (R 4.4.1)
#> rmarkdown 2.29 2024-11-04 [1] CRAN (R 4.4.1)
#> rstudioapi 0.17.1 2024-10-22 [1] CRAN (R 4.4.1)
#> scales 1.3.0 2023-11-28 [1] CRAN (R 4.4.0)
#> sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.4.0)
#> tibble 3.2.1 2023-03-20 [1] CRAN (R 4.4.0)
#> tidyr 1.3.1 2024-01-24 [1] CRAN (R 4.4.0)
#> tidyselect 1.2.1 2024-03-11 [1] CRAN (R 4.4.0)
#> utf8 1.2.4 2023-10-22 [1] CRAN (R 4.4.0)
#> vctrs 0.6.5 2023-12-01 [1] CRAN (R 4.4.0)
#> vipor 0.4.7 2023-12-18 [1] CRAN (R 4.4.0)
#> withr 3.0.2 2024-10-28 [1] CRAN (R 4.4.1)
#> xfun 0.50.5 2025-01-15 [1] Github (yihui/xfun@116d689)
#> xml2 1.3.6 2023-12-04 [1] CRAN (R 4.4.0)
#> yaml 2.3.10 2024-07-26 [1] CRAN (R 4.4.0)
#>
#> [1] /Users/nick/Library/R/arm64/4.4/library
#> [2] /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library
#>
#> ──────────────────────────────────────────────────────────────────────────────