Надежда, чтобы сделать код более эффективным R

Question

я следующий код: R

pp = function(N,J,K){ 
    for(i in 1:100){ 
     pai=runif(J) 
     alpha=matrix(rbinom(N*K,1,0.5),nrow=N) 
     Q=matrix(rbinom(J*K,1,0.5),nrow=J) 
     r=matrix(runif(J*K),nrow=J) 
     ta=r^Q 
     arrayalpha=array(rep((1-alpha), J),c(N,K,J)) 
     arrayta=array(rep(ta, N),c(J,K,N)) 
     arraytap=aperm(arrayta, c(3,2,1)) 
     tare=arraytap^arrayalpha #ta^re 
     arrayprod=apply(tare,c(1,3),prod) 
     repai=t(matrix(rep(pai,N),nrow=J)) 
     predarray=t(arrayprod*repai) 
    } 
    predarray    
} 

> system.time(pp(500,20,5)) 
    user system elapsed 
    5.381 0.008 12.468 

Как я могу сделать его более эффективным? Спасибо за помощь.

Answer 1

Учитывая, что predarray полностью заменяется на каждый цикл, вы получите те же результаты без повторения (но быстрее).

Серьезно, хотя, вы не хотите, чтобы накопить predarray ?, как:

pp = function(N,J,K){ 
    predarray <- array(numeric(100*J*N), c(100, J, N)) 
    for(i in 1:100){ 
     ... 
     predarray[i, , ] <- t(arrayprod*repai)} 
    predarray    
} 

Это было бы вернуть 100 predarray (JxN) матрицы (каждый как тусклым: 1 ломтик (100xJxN) массива). Кроме того, должно быть уменьшить ваше время с момента выделения массива (по сравнению с его распределением 100 раз, как в вашем коде).

---

Чтобы сделать его более эффективным, ваша цель должна состоять в том, чтобы избежать явного цикла R. Лучший способ сделать это - явное распараллеливание, поскольку проблема кажется «смущающей параллелью» (, т. Е. каждая итерация должна быть полностью независимой от других).

---

Другой вариант заключается в том, что, так как явно R перекручивание происходит медленнее, чем неявное C обхвата (который используется в векторизованных Операций apply семейство функций) и пакет plyr. Вы должны идти по этому, как:

• Набор num <- 100 (т.е. число «итераций», внутри определения рр)
• Сформировать каждый из (pai, alpha) как массив с дополнительным измерением длины num (т. Е. каждый элемент или срез будет эквивалентен исходным итерациям этих переменных).
• ta и остальные операции с матрицами являются узкими местами, так как вы работаете на двух матрицах каждый раз в (ta, arrayta, tare, predarray).
• Узкие места могут быть устранены путем создания (или объединения) каждой пары матричных операндов в дополнительном измерении в том же массиве (или каждый элемент pari как элемент в списке, если не одного размера), а затем маржа, представляющая «итерацию».

Пример с Q и R:

# Q and r in one array of dimensions (J x K x iterations x 2) 
# 100 Qs are stored in [,,,1] and 100 rs in [,,,2] 
Qr <- array(c(Q=rbinom(J*K*100,1,0.5), 
       r=runif(J*K*100)), 
      dim=c(J=J, K=K, iter=100, var=2)) 

# Third dimension is equivalent to each iteration 
Qr[,,1,] 

# And you operate each Q^r using apply, over each iteration 
library(package=plyr) # plyr is needed for this 
ta <- alply(.data=Qr, .margins=3, .fun=function(x) x[,,2]^x[,,1])