取R中数组的最大超维数

Question

我目前正在处理一个R中尺寸为5663x1000x100的非常大的数组。我想获得100个最大值，这将是每个5663x1000矩阵的最大值。

big_array = array(data=rnorm(566300000),dim=c(5663,1000,100))

到目前为止，我尝试过的两种方法包括for循环和application(凭直觉，这不应该是最快的方法(。

maximas = rep(0,100)
# Method 1 - Runs in 17 seconds
for(i in seq(1,100)){
maximas[i]=max(big_array[,,i])
}
# Method 2 - Runs in 36 seconds
apply(big_array,3,max)

我认为，由于数组数据结构，有一种更快的方法来运行它。我已经考虑过pmax()，但从我所看到的情况来看，我必须重塑我的数据(考虑到数组几乎是4GB，我不想创建另一个对象(。这段代码已经是正在并行化的代码的一部分了，所以我无法进一步并行化它。

任何想法都会大有帮助！

Answer 1

为什么不直接对Rcpp和RcppArmadillo执行此操作？试试这个

library(Rcpp)
library(RcppArmadillo)
cppFunction('NumericVector max_slice(const arma::cube& Q) {
int n = Q.n_slices; 
NumericVector out(n);
for (int i; i < n; i++) {
out[i] = Q.slice(i).max();
}
return out;
}', depends = "RcppArmadillo")
str(big_array)
max_slice(big_array)

输出

> str(big_array)
num [1:5663, 1:1000, 1:100] -0.282 -0.166 1.114 -0.447 -0.255 ...
> max_slice(big_array)
[1] 5.167835 4.837959 5.026354 5.211833 5.054781 5.785444 4.782578 5.169154 5.427360 5.271900 5.197460 4.994804 4.977396 5.093390 5.124796 5.221609
[17] 5.124122 4.857690 5.230277 5.217994 4.957608 5.060677 4.943275 5.382807 5.455486 5.226405 5.598238 4.942523 5.096521 5.000764 5.257607 4.843708
[33] 4.866905 5.125437 5.662431 5.224198 5.026749 5.349403 4.987372 5.228885 5.456373 5.576859 5.166118 5.124967 4.991101 5.210636 5.057471 5.005961
[49] 5.223063 5.182867 5.333683 5.528648 5.015871 4.837031 5.311825 4.981555 5.876951 5.145006 5.107017 5.252450 5.219044 5.310852 5.081958 5.210729
[65] 5.439197 5.034269 5.339251 5.567369 5.117237 5.382006 5.332199 5.032523 5.622024 5.008994 5.537377 5.279285 5.175870 5.056068 5.019422 5.616507
[81] 5.141175 4.948246 5.262170 4.961154 5.119193 4.908987 5.175458 5.328144 5.127913 5.816863 4.745966 5.507947 5.226849 5.247738 5.336941 5.134757
[97] 4.899032 5.067129 5.615639 5.118519

基准

cppFunction('NumericVector max_slice(const arma::cube& Q) {
int n = Q.n_slices; 
NumericVector out(n);
for (int i; i < n; i++) {
out[i] = Q.slice(i).max();
}
return out;
}', depends = "RcppArmadillo")
max_vapply <- function(x) vapply(seq_len(dim(x)[3]), function(i) max(x[,,i]), numeric(1))
microbenchmark::microbenchmark(
max_vapply(big_array), max_slice(big_array), 
times = 5L
)

结果

Unit: milliseconds
expr       min        lq      mean   median        uq       max neval cld
max_vapply(big_array) 4735.7055 4789.6901 5159.8319 5380.784 5428.8319 5464.1480     5   b
max_slice(big_array)  724.8582  742.0412  800.8939  747.811  833.2658  956.4935     5  a