Суммарная сумма в матрице

У меня есть матрица типа

A= [ 1 2 4
     2 3 1
     3 1 2 ]

и я бы хотел рассчитать его совокупную сумму за строкой и столбцом, то есть я хочу, чтобы результат был

B = [ 1  3  7 
      3  8  13 
      6  12 19 ]

Любые идеи о том, как быстро сделать это в R? (Возможно, используя функцию cumsum) (У меня огромные матрицы)

Спасибо!

Ответы

Ответ 1

Однострочный:

t(apply(apply(A, 2, cumsum)), 1, cumsum))

Основное наблюдение заключается в том, что вы можете сначала вычислить суммарные суммы по столбцам, а затем кумулятивную сумму этой матрицы по строкам.

Примечание. При выполнении строк вам необходимо транспонировать результирующую матрицу.

Ваш пример:

> apply(A, 2, cumsum)
     [,1] [,2] [,3]
[1,]    1    2    4
[2,]    3    5    5
[3,]    6    6    7

> t(apply(apply(A, 2, cumsum), 1, cumsum))
     [,1] [,2] [,3]
[1,]    1    3    7
[2,]    3    8   13
[3,]    6   12   19

О производительности: теперь я знаю, насколько хорошо этот подход масштабируется до больших матриц. Сложность, это должно быть близко к оптимальному. Обычно apply не так уж плох в производительности.

Изменить

Теперь мне стало любопытно, какой подход лучше? Короткий ориентир:

> A <- matrix(runif(1000*1000, 1, 500), 1000)
> 
> system.time(
+   B <- t(apply(apply(A, 2, cumsum), 1, cumsum))
+ )
       User      System     elapsed 
      0.082       0.011       0.093 
> 
> system.time(
+   C <- lower.tri(diag(nrow(A)), diag = TRUE) %*% A %*% upper.tri(diag(ncol(A)), diag = TRUE)
+ )
       User      System     elapsed 
      1.519       0.016       1.530

Таким образом: Применение превосходит матричное умножение в 15 раз. (Только для сравнения: MATLAB необходимо 0.10719 секунд.) Результаты не удивляют, так как apply -версия может быть выполнена в O (n ^ 2), в то время как для матричного умножения потребуется ок. O (n ^ 2.7). Таким образом, все оптимизации, предлагаемые предложением матриц, должны быть потеряны, если n достаточно велико.

Ответ 2

Вот более эффективная реализация с использованием пакета matrixStats и более крупной примерной матрицы:

library(matrixStats)
A <- matrix(runif(10000*10000, 1, 500), 10000)

# Thilo answer
system.time(B <- t(apply(apply(A, 2, cumsum), 1, cumsum)))
user  system elapsed 
3.684   0.504   4.201

# using matrixStats
system.time(C <- colCumsums(rowCumsums(A)))
user  system elapsed 
0.164   0.068   0.233 

all.equal(B, C)
[1] TRUE