2014-02-27 81 views
0

我試圖從我經常使用的包加速R函數,所以任何幫助矢量化下面的for循環將非常感激!需要幫助矢量化for循環在R

y <- array(0, dim=c(75, 12)) 
samp <- function(x) x<-sample(c(0,1), 1) 
y <- apply(y, c(1,2), samp) 

nr <- nrow(y) 
nc <- ncol(y) 
rs <- rowSums(y) 
p <- colSums(y) 
out <- matrix(0, nrow = nr, ncol = nc) 

for (i in 1:nr) { 
    out[i, sample.int(nc, rs[i], prob = p)] <- 1 
} 

我遇到困難的問題是循環內對象'rs'的引用。

有什麼建議嗎?

+1

有一個[sample'樣式的RcppArmadillo實現](http://gallery.rcpp.org/articles/using-the-Rcpp-based-sample-implementation/)。所以,你可以嘗試用Rcpp來實現,看看它是否更快。 – Roland

回答

1

這裏有兩種選擇:

這其中使用了有點氣餒<<-操作:

lapply(1:nr, function(i) out[i, sample.int(nc, rs[i], prob = p)] <<- 1) 

這一次使用較爲傳統的索引:

out[do.call('rbind',sapply(1:nr, function(i) cbind(i,sample.int(nc, rs[i], prob = p))))] <- 1 

我想你也可以使用Vectorize在你的功能上做一個隱含的mapply

z <- Vectorize(sample.int, vectorize.args='size')(nc, rs, prob=p) 
out[cbind(rep(1:length(z), sapply(z, length)), unlist(z))] <- 1 

但我不認爲這一定是更清潔。

而且,事實上,@Roland是正確的,所有的這些都不僅僅是做for循環較慢:

> microbenchmark(op(), t1(), t2(), t3()) 
Unit: microseconds 
expr  min  lq median  uq  max neval 
op() 494.970 513.8290 521.7195 532.3040 1902.898 100 
t1() 591.962 602.1615 609.4745 617.5570 2369.385 100 
t2() 734.756 754.7700 764.3925 782.4825 2205.421 100 
t3() 642.383 672.9815 711.4700 763.8150 2283.169 100 

耶自由利益混淆!

+1

使用'lapply'和'<< - '來代替'for'應該不會受到某種程度的阻礙。 –

+0

我不認爲你的建議比原來的'for'循環更快。他們只是更混亂。 – Roland

+0

@羅蘭我同意這一點。 – Thomas