2014-02-27 55 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