R:向量中不同维数的矩阵元素之和

Question

我有一个包含许多矩阵的列表，这些矩阵对应于计数分布中的重复绘制(列表示原始计数数据，行表示模拟计数数据，单元格值是来自一次模拟的观察计数)。

我希望能够将列表中的矩阵求和，形成一个"超级"。矩阵，包含模拟中观测值的计数。我知道矩阵不应该相加，除非它们是相同的维度，所以如果有人有一个想法如何使用像grid.expand这样的东西，所有的矩阵都是相同的大小，我愿意接受。

我还应该注意到，我并不局限于矩阵——我也可以提取大小不等的计数向量，但我仍然不知道如何(在一个矢量化的解决方案中)对许多大小不等的计数向量求和而不循环。这可能已经在某个地方得到了回答，但乍一看，我找不到任何有意义的东西在许多模拟的矢量化解决方案中。

我在下面做了一个玩具表示-mtcars$carb可能不是理想的变量，但我想要一个小的泊松分布变量，所以请坚持我。

感谢任何人提供的任何解决方案!

library(tidyverse)
draw <- mtcars %>% 
select(carb) %>% 
mutate(carb2 = NA)
# initialize storage vector to hold results from 5 simulations
res <- vector(mode = "list", length = 5)
# run simulation 5 times
for (i in 1:5) { 

# draw new counts once from carb ~ Poisson in carb2 with lambda = mean(carb)
for (j in 1:nrow(draw)) {
draw$carb2[j] <- rpois(1, mean(draw$carb))
}

# get "crosstabs" matrix of obs counts 
# cols are original values from carb, rows are sim values from carb2; 
# values are counts
sim <- draw %>% 
group_by(carb, carb2) %>% 
tally() %>% 
spread(carb, n) %>% 
replace(is.na(.), 0) %>%   # replace NA with 0
as.matrix()

res[[i]] <- sim

}

输出如下:

> res
[[1]]
carb2 1 2 3 4 6 8
[1,]     0 1 1 0 2 0 0
[2,]     1 1 0 0 3 0 0
[3,]     2 1 4 1 2 0 1
[4,]     3 3 1 0 0 1 0
[5,]     4 0 3 0 2 0 0
[6,]     5 1 1 2 1 0 0
[[2]]
carb2 1 2 3 4 6 8
[1,]     0 0 2 0 1 0 0
[2,]     1 1 0 0 1 1 0
[3,]     2 3 1 1 1 0 0
[4,]     3 2 2 0 3 0 1
[5,]     4 0 3 1 3 0 0
[6,]     5 0 1 1 1 0 0
[7,]     6 0 1 0 0 0 0
[8,]     7 1 0 0 0 0 0
[[3]]
carb2 1 2 3 4 6 8
[1,]     1 1 5 0 2 0 0
[2,]     2 2 2 1 3 0 0
[3,]     3 3 2 1 1 0 1
[4,]     4 1 1 0 2 1 0
[5,]     5 0 0 0 1 0 0
[6,]     6 0 0 1 0 0 0
[7,]     8 0 0 0 1 0 0
[[4]]
carb2 1 2 3 4 6 8
[1,]     0 1 0 0 1 0 0
[2,]     1 0 1 0 1 0 0
[3,]     2 2 2 1 1 0 0
[4,]     3 0 2 0 5 1 1
[5,]     4 1 2 2 0 0 0
[6,]     5 3 3 0 2 0 0
[[5]]
carb2 1 2 3 4 6 8
[1,]     0 1 0 0 0 0 0
[2,]     1 1 3 1 1 0 0
[3,]     2 2 3 1 2 1 0
[4,]     3 1 0 0 3 0 1
[5,]     4 1 3 0 3 0 0
[6,]     5 1 0 1 0 0 0
[7,]     6 0 1 0 1 0 0

我想把上面五个模拟的矩阵求和成下面的格式:

carb2   1   2   3   4   6   8
0   2   1   0   2   0   0
1   5   11  1   10  1   0
2   10  12  5   9   1   1
3   9   7   1   12  2   4
4   3   12  3   10  1   0
5   5   5   4   5   0   0
6   0   2   1   1   0   0
7   1   0   0   0   0   0
8   0   0   0   1   0   0

Answer 1

有个办法。

主要思想是将res的数据结构从列表改为3d数组。前两个维度分别为carb2和carb，第三个维度为模拟次数。事先没有办法猜测rpois的上限，所以在下面的代码中，我选择了0.999分位数，希望它足够了。然后将每次迭代的模拟结果sim赋值给行与carb2匹配的数组切片。最后，添加所有的切片。

suppressPackageStartupMessages(
library(tidyverse)
)
set.seed(2022)
draw <- mtcars %>% 
select(carb) %>% 
mutate(carb2 = NA)
cols <- c("carb2", sort(unique(mtcars$carb)))
rows <- 0:qpois(0.999, mean(draw$carb))
R <- 5L
res <- array(NA_integer_, dim = c(length(rows), length(cols), R),
dimnames = list(rows, cols, 1:R))
for(i in 1:R) {
draw$carb2 <- rpois(nrow(draw), mean(draw$carb))
sim <- draw %>% 
group_by(carb, carb2) %>% 
tally() %>% 
spread(carb, n) %>% 
replace(is.na(.), 0) %>%   # replace NA with 0
as.matrix()
j <- match(sim[, "carb2"], rows)
res[j, , i] <- sim
}
cbind(
carb2 = seq_len(nrow(res)) - 1L,
apply(res[, -1, ], 1:2, sum, na.rm = TRUE)
)
#>   carb2  1  2 3  4 6 8
#> 0     0  0  3 0  3 1 0
#> 1     1  6 11 5  6 0 0
#> 2     2  7  9 1  9 1 1
#> 3     3 11 12 4 13 1 2
#> 4     4  5  6 2 10 2 2
#> 5     5  4  6 2  7 0 0
#> 6     6  2  1 0  1 0 0
#> 7     7  0  1 1  1 0 0
#> 8     8  0  0 0  0 0 0
#> 9     9  0  1 0  0 0 0

^{由reprex包(v2.0.1)创建于2022-07-14}