在R中对大矩阵的每一行进行排序的最快方法

Question

我有一个大矩阵：

set.seed(1)
a <- matrix(runif(9e+07),ncol=300)

我想对矩阵中的每一行进行排序：

> system.time(sorted <- t(apply(a,1,sort)))
   user  system elapsed 
  42.48    3.40   45.88 

我有很多 RAM 可以使用，但我想要一种更快的方法来执行此操作。

Answer 1

好吧，我不知道在 R 中排序更快的方法有多少，问题是您只对 300 个值进行排序，但是很多次。不过，您可以通过直接调用sort.int 并使用method='quick' 来获得一些额外的性能：

set.seed(1)
a <- matrix(runif(9e+07),ncol=300)

# Your original code
system.time(sorted <- t(apply(a,1,sort))) # 31 secs

# sort.int with method='quick'
system.time(sorted2 <- t(apply(a,1,sort.int, method='quick'))) # 27 secs

# using a for-loop is slightly faster than apply (and avoids transpose):
system.time({sorted3 <- a; for(i in seq_len(nrow(a))) sorted3[i,] <- sort.int(a[i,], method='quick') }) # 26 secs

但更好的方法应该是使用并行包对矩阵的各个部分进行并行排序。但是，传输数据的开销似乎太大了，在我的机器上它开始交换，因为我“只有”8 GB 内存：

library(parallel)
cl <- makeCluster(4)
system.time(sorted4 <- t(parApply(cl,a,1,sort.int, method='quick'))) # Forever...
stopCluster(cl)

Answer 2

grr 包包含一个替代排序方法，可用于加速此特定操作（我已稍微减小了矩阵大小，以便此基准测试不会永远持续）：

> set.seed(1)
> a <- matrix(runif(9e+06),ncol=300)
> microbenchmark::microbenchmark(sorted <- t(apply(a,1,sort))
+                                ,sorted2 <- t(apply(a,1,sort.int, method='quick'))
+                                ,sorted3 <- t(apply(a,1,grr::sort2)),times=3,unit='s')
Unit: seconds
                                                  expr       min       lq     mean   median       uq      max neval
                        sorted <- t(apply(a, 1, sort)) 1.7699799 1.865829 1.961853 1.961678 2.057790 2.153902     3
 sorted2 <- t(apply(a, 1, sort.int, method = "quick")) 1.6162934 1.619922 1.694914 1.623551 1.734224 1.844898     3
                 sorted3 <- t(apply(a, 1, grr::sort2)) 0.9316073 1.003978 1.050569 1.076348 1.110049 1.143750     3

当矩阵包含字符时，差异变得显着：

> set.seed(1)
> a <- matrix(sample(letters,size = 9e6,replace = TRUE),ncol=300)
> microbenchmark::microbenchmark(sorted <- t(apply(a,1,sort))
+                                ,sorted2 <- t(apply(a,1,sort.int, method='quick'))
+                                ,sorted3 <- t(apply(a,1,grr::sort2)),times=3)
Unit: seconds
                                                  expr       min        lq      mean    median        uq      max neval
                        sorted <- t(apply(a, 1, sort)) 15.436045 15.479742 15.552009 15.523440 15.609991 15.69654     3
 sorted2 <- t(apply(a, 1, sort.int, method = "quick")) 15.099618 15.340577 15.447823 15.581536 15.621925 15.66231     3
                 sorted3 <- t(apply(a, 1, grr::sort2))  1.728663  1.733756  1.780737  1.738848  1.806774  1.87470     3

三个结果都是一样的。

> identical(sorted,sorted2,sorted3)
[1] TRUE

Answer 3

Martin Morgan 的另一个出色方法，在 Fastest way to select i-th highest value from row and assign to new column 中没有使用任何外部包：

matrix(a[order(row(a), a)], ncol=ncol(a), byrow=TRUE)

在同一个链接的 cmets 下还有一个等价的按列排序。

使用与 Craig 相同的数据的计时代码：

set.seed(1)
a <- matrix(runif(9e7),ncol=300)

use_for <- function(){
    sorted3 <- a
    for(i in seq_len(nrow(a))) 
        sorted3[i,] <- sort.int(a[i,], method='quick') 
    sorted3
}

microbenchmark::microbenchmark(times=3L,
    t(apply(a,1,sort)),
    t(apply(a,1,sort.int, method='quick')),
    use_for(),
    Rfast::rowSort(a),
    t(apply(a,1,grr::sort2)),
    mmtd=matrix(a[order(row(a), a)], ncol=ncol(a), byrow=TRUE)
)

时间安排：

Unit: seconds
                                       expr       min        lq      mean    median        uq       max neval
                       t(apply(a, 1, sort)) 24.233418 24.305339 24.389650 24.377260 24.467766 24.558272     3
 t(apply(a, 1, sort.int, method = "quick")) 17.024010 17.156722 17.524487 17.289433 17.774726 18.260019     3
                                  use_for() 13.384958 13.873367 14.131813 14.361776 14.505241 14.648705     3
                          Rfast::rowSort(a)  3.758765  4.607609  5.136865  5.456452  5.825914  6.195377     3
                 t(apply(a, 1, grr::sort2))  9.810774  9.955199 10.310328 10.099624 10.560106 11.020587     3
                                       mmtd  6.147010  6.177769  6.302549  6.208528  6.380318  6.552108     3

---

为了呈现更完整的画面，对字符类进行另一次测试（不包括Rfast::rowSort，因为它无法处理字符类）：

set.seed(1)
a <- matrix(sample(letters, 9e6, TRUE),ncol=300)

microbenchmark::microbenchmark(times=1L,
    t(apply(a,1,sort)),
    t(apply(a,1,sort.int, method='quick')),
    use_for(),
    #Rfast::rowSort(a),
    t(apply(a,1,grr::sort2)),
    mmtd=matrix(a[order(row(a), a, method="radix")], ncol=ncol(a), byrow=TRUE)
)

时间安排：

Unit: milliseconds
                                       expr        min         lq       mean     median         uq        max neval
                       t(apply(a, 1, sort)) 14848.4356 14848.4356 14848.4356 14848.4356 14848.4356 14848.4356     1
 t(apply(a, 1, sort.int, method = "quick")) 15061.0993 15061.0993 15061.0993 15061.0993 15061.0993 15061.0993     1
                                  use_for() 14144.1264 14144.1264 14144.1264 14144.1264 14144.1264 14144.1264     1
                 t(apply(a, 1, grr::sort2))  1831.1429  1831.1429  1831.1429  1831.1429  1831.1429  1831.1429     1
                                       mmtd   440.9158   440.9158   440.9158   440.9158   440.9158   440.9158     1

---

正面交锋：

set.seed(1)
a <- matrix(sample(letters, 9e7, TRUE),ncol=300)
microbenchmark::microbenchmark(times=1L,
    t(apply(a,1,grr::sort2)),
    mmtd=matrix(a[order(row(a), a, method="radix")], ncol=ncol(a), byrow=TRUE)
)

时间安排：

Unit: seconds
                       expr       min        lq      mean    median        uq       max neval
 t(apply(a, 1, grr::sort2)) 19.273225 19.273225 19.273225 19.273225 19.273225 19.273225     1
                       mmtd  3.854117  3.854117  3.854117  3.854117  3.854117  3.854117     1

R 版本：

R version 4.0.3 (2020-10-10)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 18363)