如何在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