-2

我必须在一个大矩阵上计算一系列统计数据，并且我想以最有效的方式使用向量作为分组因子来做到这一点。

行是我要分组的变量，而列是样本。

例如：

mat = matrix(seq(1,10000), ncol  = 100)
vect_group = c(1,1,1,1,1,2,2,2,3,3,3, ...)

我想计算索引为 1、2、3 等的所有行的列的平均值。因此，在这种情况下，获得一个新矩阵，其行vect_group数与匹配列中的级别和相应统计量一样多。

到目前为止，我通过索引获得了这个循环，并且每次都在这些子矩阵上使用 apply，但我想加快这个过程。我尝试过doParallel，foreach但没有成功。

我正在努力解决的关键部分是生成较小矩阵的拆分/聚合过程。另外，我不知道开销是否会影响多线程计算的选择。

4

2 回答 2

1

我不知道你是否需要多线程。

我测试了两种解决方案，一种使用 base R，另一种使用dplyr. 两者在基准测试中都非常快。

mat <- matrix(seq(1,10000), ncol  = 100)
vect_group <- rep(1:10, each = 10)

#--
library(dplyr)

#-- Base R
splitData <- split(as.data.frame(mat), vect_group)
meansPerGroup <- sapply(splitData, colMeans)

#-- Dplyr
df <- data.frame(mat, vect_group)
meansPerGroup <- df %>%
    group_by(vect_group) %>%
    summarize_at(vars(colnames(mat)), mean)

然后我对这两种解决方案进行了基准测试：

rbenchmark::benchmark(replications = 5000,
    baseR = function(mat = mat, vect_group = vect_group) {
        splitData <- split(as.data.frame(mat), vect_group)
        meansPerGroup <- sapply(splitData, colMeans)
    },
    dplyr = function(df = df, vect_group = vect_group) {
        meansPerGroup <- df %>%
            group_by(vect_group) %>%
            summarize_at(vars(colnames(mat)), mean)
    })

基准测试结果：

   test replications elapsed relative user.self sys.self user.child sys.child
1 baseR         5000   0.006      1.2     0.006        0          0         0
2 dplyr         5000   0.005      1.0     0.006        0          0         0
于 2019-02-07T16:28:32.437 回答
0

我同意@csgroen 的观点，即并行执行此计算可能是不必要的，因为计算平均值非常快并且设置它需要开销，但这可能取决于您问题的规模。你的矩阵有多大？

可能不并行的最快方法是使用data.table. 我在下面对一些方法进行了基准测试，包括上一个答案（尽管我无法让 dplyr 版本在我的计算机上运行——我认为是因为mat没有列名）。Data.table 平均需要大约 3 毫秒，并且聚合也不远了。

mat <-  matrix(seq(1,10000), ncol  = 100)
vect_group  = rep(1:10, each = 10)

fn1_agg <- function(mat, vg) {
  aggregate(c(mat)~rep(vg, ncol(mat)), FUN = mean)
}

fn2_dt <- function(mat, vg){
  DT <- data.table::data.table(m = c(mat), v = rep(vg, ncol(mat)))
  data.table::setkey(DT, v)
  DT[, list(m = mean(m)), by = v]
}

fn3_split <- function(mat, vg) {
  splitData <- split(as.data.frame(mat), vect_group)
  sapply(splitData, colMeans)
}

microbenchmark::microbenchmark(fn1_agg(mat, vect_group),
                               fn2_dt(mat, vect_group),
                               fn3_split(mat, vect_group))
#> Unit: milliseconds
#>                        expr       min        lq      mean    median
#>    fn1_agg(mat, vect_group)  5.169709  5.437589  6.122462  6.293567
#>     fn2_dt(mat, vect_group)  1.197218  1.291972  3.004166  1.472097
#>  fn3_split(mat, vect_group) 15.480264 15.751230 16.998514 16.267098
#>         uq        max neval cld
#>   6.481626   9.454458   100  b 
#>   1.538948 142.368800   100 a  
#>  17.060969  60.686907   100   c

由reprex 包（v0.2.1）于 2019 年 2 月 7 日创建

于 2019-02-07T16:49:35.253 回答