3

发布在 R、NULL 和 NA 中分配矩阵的最佳方法？表明在 R 中编写自己的矩阵分配函数可以比使用 R 的内置 matrix() 函数预分配大矩阵快 8 到 10 倍。

有谁知道为什么手工制作的功能要快得多？R 在如此慢的 matrix() 内部做了什么？谢谢。

这是我系统上的代码：

create.matrix <- function( nrow, ncol ) {
x<-matrix()
length(x) <- nrow*ncol
dim(x) <- c(nrow,ncol)
x
}

system.time( x <- matrix(nrow=10000, ncol=9999) )
user  system elapsed 
1.989   0.136   2.127 

system.time( y <- create.matrix( 10000, 9999 ) )
user  system elapsed 
0.192   0.141   0.332 
identical(x,y)
[1] TRUE

我向那些评论认为用户定义的函数速度较慢的人道歉，因为上面链接的答案中发布的内容不一致。我正在查看用户时间，在上面的链接中大约快 8 倍，而在我的系统上，用户定义的与内置的大约快 10 倍。

回应 Joshua 的会话信息请求：

> sessionInfo()
R version 2.12.1 (2010-12-16)
Platform: x86_64-unknown-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=C              LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

loaded via a namespace (and not attached):
[1] tools_2.12.1

此外，我尝试运行 Simon 的三个示例，而 Simon 给出的第三个示例最快，对我来说却是最慢的：

> system.time(matrix(NA, nrow=10000, ncol=9999)) 
   user  system elapsed 
  2.011   0.159   2.171 
> system.time({x=NA; length(x)=99990000; dim(x)=c(10000,9999); x}) 
   user  system elapsed 
  0.194   0.137   0.330 
> system.time(matrix(logical(0), nrow=10000, ncol=9999)) 
   user  system elapsed 
  4.180   0.200   4.385

然而，我仍然认为 Simon 的想法可能是正确的，即matrix()最初分配一个 1x1 矩阵然后复制它。有人知道关于 R 内部的任何好的文档吗？谢谢。

4

4 回答 4

8

问题是你的matrix电话比你想象的要复杂一些。比较以下版本：

# copy NA matrix
> system.time(matrix(NA, nrow=10000, ncol=9999))
   user  system elapsed 
  1.272   0.224   1.496 

# replicate NA vector (faster version of what you used)
> system.time({x=NA; length(x)=99990000; dim(x)=c(10000,9999); x})
   user  system elapsed 
  0.292   0.260   0.552 

# fastest - just allocate a matrix filled with NAs 
> system.time(matrix(logical(0), nrow=10000, ncol=9999))
   user  system elapsed 
  0.184   0.308   0.495

因此，在您的示例中，您实际上是在创建一个 1 x 1 NA 矩阵，该矩阵被复制到您指定的大小 - 最慢的方法。对向量执行相同操作会更快（因为它不需要对列使用模）-您以一种有点复杂的方式进行操作（通过创建矩阵，将其转换为向量，然后再转换回矩阵），但是这个想法是一样的。最后，如果你只使用一个空向量，那么矩阵将简单地用NA你想要的东西填充，因此不需要额外的工作（最快）。

编辑一个重要的说明，虽然：马修的建议是正确的，虽然它没有涉及（因为他引用的代码是logical(0)案例，而不是NA案例）。无意中我在上述时间运行了 R-devel，所以发布的 R 中的时间会有所不同。

于 2012-09-03T03:06:52.683 回答
5

我将对这些评论提出异议，尽管我确实理解其中的大部分内容。问题在于，引用的帖子的答案存在内部矛盾，评论者一直依赖而没有检查。用户和系统的时间没有按应有的方式正确加起来。

 create.matrix <- function(size) {
  x <- matrix()
  length(x) <- size^2
  dim(x) <- c(size,size)
  x
  }
  system.time(x <- matrix(data=NA,nrow=10000,ncol=10000))
#   user  system elapsed 
#  0.464   0.226   0.688 
 system.time(y <- create.matrix(size=10000))
#   user  system elapsed 
#  0.177   0.239   0.414

我怀疑效率实际上是通过用户定义的函数只能创建方阵和“矩阵”需要对更一般情况下参数的有效性进行检查这一事实来实现的。

编辑：我看到你已经反驳了我的一个假设（关于方阵限制），我还会注意到我的另一个假设，即这是由于惰性评估导致的，我的测试也失败了。这种差异确实没有意义，因为用户代码使用了该matrix功能。

于 2012-08-31T19:00:14.787 回答
4

不确定这是不是原因（可能是不同的低效率），但do_matrix在 src/array.c 中有一个类型开关，其中包含：

case LGLSXP :
    for (i = 0; i < nr; i++)
    for (j = 0; j < nc; j++)
        LOGICAL(ans)[i + j * NR] = NA_LOGICAL;

这看起来是页面效率低下。认为应该是：

case LGLSXP :
    for (j = 0; j < nc; j++)
    for (i = 0; i < nr; i++)
        LOGICAL(ans)[i + j * NR] = NA_LOGICAL;

或更简单地说：

case LGLSXP :
    for (i = 0; i < nc*nr; i++)
        LOGICAL(ans)[i] = NA_LOGICAL;

（需要一些微调，因为is NRtype R_xlen_twhilei和are type ）。ncnrint

更新：

发布到 r-devel 后：

array.c 中的 do_matrix 中可能的页面效率低下

Simon Urbanek 现在已经承诺对 R 进行更改。它现在使用上面的单一索引方法：

array.c 的最新实时版本

但正如西蒙所说，我在上面介绍了自己，这似乎并不能解决问题提出的特定问题。第二个，不同的，低效率也需要被发现和修复。

这可能是后续修复的内容。这结合了新代码的页面效率（现在在 R 中），但是，切换到使用时matrix(data=NA)（R 的默认值）。通过避免这种情况，这避免了copyMatrix西蒙在他的回答中提到的模数。copyMatrixNA

目前，array.c 中的 do_matrix 具有：

if(lendat) {
    if (isVector(vals))
        copyMatrix(ans, vals, byrow);
    else
        copyListMatrix(ans, vals, byrow);
} else if (isVector(vals)) { 
    // fill with NAs in the new page efficient way that Simon already committed.
}

可能如下。我不知道ISNA()C 级别的函数需要SEXP输入，所以已经编写了那个长手（西蒙，有更好的方法吗？）：

if(lendat && // but not a single NA, basically :
             !(lendat==1 &&
                  ((isLogical(vals) && LOGICAL(vals)[0] == NA_LOGICAL) ||
                   (isReal(vals) && ISNA(REAL(vals)[0])) ||
                   (isInteger(vals) && INTEGER(vals)[0] == NA_INTEGER)))) {
    if (isVector(vals))
        copyMatrix(ans, vals, byrow);
    else
        copyListMatrix(ans, vals, byrow);
} else if (isVector(vals)) { 
    // fill with NAs in the new page efficient way that Simon already committed.
    // this branch will now run when dat is a single NA, too
}
于 2012-09-03T01:11:09.510 回答
1

唔。是的，这很奇怪。...而且这仍然稍微快一点 - 尽管它更像 matrix() ，因为它允许单个数据参数（但它必须是标量）：

create.matrix2 <- function(data=NA, nrow, ncol) {
  x <- rep.int(data[[1]], nrow*ncol)
  dim(x) <- c(nrow, ncol)
  x
}
system.time( x <- matrix(nrow=10000, ncol=9999) ) # 0.387 secs
system.time( y <- create.matrix(nrow=10000, ncol=9999) )  # 0.199 secs
system.time( z <- create.matrix2(nrow=10000, ncol=9999) ) # 0.173 secs
identical(x,z) # TRUE

...我猜用于创建矩阵的内部代码会浪费一些东西（或者可能有用，但我想不出那会是什么）...

哦，因为它处理data任何长度，它最终可能会得到类似于rep(data, length.out=nrow*ncol)相当慢的东西：

system.time( rep(NA, length.out=10000*9999) ) # 1.5 secs!

无论如何，绝对有改进的余地！

于 2012-08-31T21:17:54.223 回答