Суммирование строк матрицы на основе индекса столбца

Я пытаюсь перейти от матрицы, в которой есть столбцы, которые «принадлежат друг другу», к той, в которой были сформированы суммы строк соответствующих подматриц. т.е. идущий от

     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16]
[1,]    1    5    9   13   17   21   25   29   33    37    41    45    49    53    57    61
[2,]    2    6   10   14   18   22   26   30   34    38    42    46    50    54    58    62
[3,]    3    7   11   15   19   23   27   31   35    39    43    47    51    55    59    63
[4,]    4    8   12   16   20   24   28   32   36    40    44    48    52    56    60    64

to

     [,1] [,2] [,3] [,4] [,5]
[1,]   15   30   46  185  220
[2,]   18   32   48  190  224
[3,]   21   34   50  195  228
[4,]   24   36   52  200  232

Я предполагаю, что должен быть какой-то гораздо более элегантный и быстрый способ сделать это, чем перебирать индексы, как я делаю ниже (в частности, моя реальная матрица будет больше похожа на 4000 на многие тысячи).

example <- matrix(1:64, nrow=4) myindex <- c(1,1,1,2,2,3,3,4,4,4,4,4,5,5,5,5) summed <- matrix( rep(unique(myindex), each=dim(example)[1]), nrow=dim(example)[1]) for (i in 1:length(unique(myindex))){ summed[,i] <- apply(X=example[,(myindex==i)], MARGIN=1, FUN=sum) }

Вероятно, это отсутствие у меня опыта работы с apply и tapply мешает мне понять это. Быстрый подход dplyr, конечно, также приветствуется.


person Björn    schedule 28.03.2018    source источник


Ответы (4)


Мы можем использовать один лайнер с sapply:

sapply(unique(myindex), function(x) rowSums(example[, which(myindex == x), drop = FALSE]))

     [,1] [,2] [,3] [,4] [,5]
[1,]   15   30   46  185  220
[2,]   18   32   48  190  224
[3,]   21   34   50  195  228
[4,]   24   36   52  200  232

Мы позволяем sapply перебирать все уникальные значения myindex и используем which для определения столбцов, которые должны быть включены в rowSums.


Изменить: включено drop = FALSE, чтобы предотвратить упрощение отдельных индексов до вектора. Спасибо @mt1022 за указание на ошибку!

person LAP    schedule 28.03.2018
comment
Спасибо за один из быстрых ответов на мой пример 4000 на 3020, он кажется самым быстрым из 3 предложенных ответов, поэтому я отмечу его как принятый ответ. - person Björn; 28.03.2018

Мы также можем сделать это, splitting

sapply(split.default(as.data.frame(example), myindex), rowSums)
#     1  2  3   4   5
#[1,] 15 30 46 185 220
#[2,] 18 32 48 190 224
#[3,] 21 34 50 195 228
#[4,] 24 36 52 200 232
person akrun    schedule 28.03.2018

Другой подход...

example <- matrix(1:64, nrow=4)
myindex <- c(1,1,1,2,2,3,3,4,4,4,4,4,5,5,5,5)

summed <- t(apply(example,1,cumsum))
summed <- summed[,cumsum(rle(myindex)$lengths)]
summed[,-1] <- t(apply(summed,1,diff))
summed

     [,1] [,2] [,3] [,4] [,5]
[1,]   15   30   46  185  220
[2,]   18   32   48  190  224
[3,]   21   34   50  195  228
[4,]   24   36   52  200  232
person Andrew Gustar    schedule 28.03.2018

Альтернативный подход с матричным умножением (менее эффективный для большого набора данных):

x <- matrix(0, nrow = ncol(example), ncol = max(myindex))
x[cbind(1:ncol(example), myindex)] <- 1
example %*% x

#      [,1] [,2] [,3] [,4] [,5]
# [1,]   15   30   46  185  220
# [2,]   18   32   48  190  224
# [3,]   21   34   50  195  228
# [4,]   24   36   52  200  232

Вот эталонный тест с примерными данными, соответствующими реальному размеру данных:

library(microbenchmark)

n_row <- 4000
n_col <- 3020
example <- matrix(rnorm(n_row * n_col), nrow = n_row)
myindex <- ceiling((1:n_col)/5)

microbenchmark(
    matrix = {
        x <- matrix(0, nrow = ncol(example), ncol = max(myindex))
        x[cbind(1:ncol(example), myindex)] <- 1
        example %*% x
    },
    split = {  # by akrun
        sapply(split.default(as.data.frame(example), myindex), rowSums)
    },
    which = {  # by LAP
        sapply(unique(myindex), function(x) rowSums(example[, which(myindex == x)]))
    },
    times = 10
)

# Unit: milliseconds
#    expr       min        lq     mean    median       uq      max neval
#  matrix 982.55727 989.65177 992.7295 992.91230 997.3704 999.0066    10
#   split 162.13377 162.57711 194.5668 167.92963 182.5335 403.8740    10
#   which  90.28227  94.82681 119.3977  96.03701 103.1125 316.9170    10
person mt1022    schedule 28.03.2018
comment
Спасибо. Пакет microbench действительно интересен. Interestingly with my real example (4000 by 3020), I ended up with Unit: milliseconds expr lq mean median uq matrix 13040.5662 13535.0503 13874.2902 14135.7251 split 2379.6066 3173.5876 2631.1031 3371.2384 which 204.2357 322.9715 254.2363 383.8717 - person Björn; 28.03.2018
comment
@Björn, кажется, что матричный подход слишком медленный для реального набора данных. Я думал, что это должно быть быстрее. - person mt1022; 28.03.2018