如何在R中更快地运行循环？

Question

我想使用arms()每次获取一个样本，并在我的函数中创建一个如下所示的循环。它运行非常缓慢。我怎么能让它跑得更快？谢谢。

library(HI)  
dmat <- matrix(0, nrow=100,ncol=30) 
system.time(
    for (d in 1:100){ 
     for (j in 1:30){ 
      y <- rep(0, 101) 
      for (i in 2:100){ 

       y[i] <- arms(0.3, function(x) (3.5+0.000001*d*j*y[i-1])*log(x)-x, 
        function(x) (x>1e-4)*(x<20), 1)  
      } 
     dmat[d, j] <- sum(y) 
     } 
    } 
) 

Answer 1

这是一个版本的基础上汤米的答案，但避免了所有的循环：

library(multicore) # or library(parallel) in 2.14.x 
set.seed(42) 
m = 100 
n = 30 
system.time({ 
    arms.C <- getNativeSymbolInfo("arms")$address 
    bounds <- 0.3 + convex.bounds(0.3, dir = 1, function(x) (x>1e-4)*(x<20)) 
    if (diff(bounds) < 1e-07) stop("pointless!") 
    # create the vector of z values 
    zval <- 0.00001 * rep(seq.int(n), m) * rep(seq.int(m), each = n) 
    # apply the inner function to each grid point and return the matrix 
    dmat <- matrix(unlist(mclapply(zval, function(z) 
      sum(unlist(lapply(seq.int(100), function(i) 
       .Call(arms.C, bounds, function(x) (3.5 + z * i) * log(x) - x, 
         0.3, 1L, parent.frame()) 
      ))) 
     )), m, byrow=TRUE) 
}) 

在多核机器，这将是非常快，因为它在核传播的负载。在单核机器上（或者对于Windows用户较差的用户），您可以使用lapply替换上面的mclapply，与Tommy的答案相比，只能稍微提高速度。但请注意，并行版本的结果会有所不同，因为它将使用不同的RNG序列。

请注意，任何需要评估R函数的C代码本质上都很慢（因为解释代码很慢）。我添加了arms.C只是为了消除所有R-> C开销以使moli高兴;）但它没有任何区别。

通过使用列专业处理（问题代码为row-major，需要重新复制，因为R矩阵始终为列专业），您可以挤出更多的毫秒。

编辑：我注意到，摩力变化不大，因为汤米回答了这个问题 - 这样，而不是你必须使用一个循环，因为y[i]都依赖于sum(...)部分，所以function(z)看起来像

function(z) { y <- 0 
    for (i in seq.int(99)) 
     y <- y + .Call(arms.C, bounds, function(x) (3.5 + z * y) * log(x) - x, 
         0.3, 1L, parent.frame()) 
    y } 

Answer 2

为什么不喜欢这个？

dat <- expand.grid(d=1:10, j=1:3, i=1:10) 

arms.func <- function(vec) { 
    require(HI) 
    dji <- vec[1]*vec[2]*vec[3] 
    arms.out <- arms(0.3, 
        function(x,params) (3.5 + 0.00001*params)*log(x) - x, 
        function(x,params) (x>1e-4)*(x<20), 
        n.sample=1, 
        params=dji) 

    return(arms.out) 
} 

dat$arms <- apply(dat,1,arms.func) 

library(plyr) 
out <- ddply(dat,.(d,j),summarise, arms=sum(arms)) 

matrix(out$arms,nrow=length(unique(out$d)),ncol=length(unique(out$j))) 

但是，它仍然单核心和耗时。但这不是R慢，它的武器功能。

Answer 3

那么，一个有效的方法是摆脱arms内部的开销。它会进行一些检查并每次调用indFunc，即使结果总是与您的情况相同。 其他一些评估也可以在循环之外完成。这些优化将我的机器上的时间从54秒减少到约6.3秒。 ......答案是一样的。

set.seed(42) 
#dmat2 <- ##RUN ORIGINAL CODE HERE## 

# Now try this: 
set.seed(42) 
dmat <- matrix(0, nrow=100,ncol=30) 
system.time({ 
    e <- new.env() 
    bounds <- 0.3 + convex.bounds(0.3, dir = 1, function(x) (x>1e-4)*(x<20)) 
    f <- function(x) (3.5+z*i)*log(x)-x 
    if (diff(bounds) < 1e-07) stop("pointless!") 
    for (d in seq_len(nrow(dmat))) { 
     for (j in seq_len(ncol(dmat))) { 
      y <- 0 
      z <- 0.00001*d*j 
      for (i in 1:100) { 
       y <- y + .Call("arms", bounds, f, 0.3, 1L, e) 
      } 
      dmat[d, j] <- y 
     } 
    } 
}) 

all.equal(dmat, dmat2) # TRUE