加快R中的循环操作
我在R中有一个很大的性能问题。我写了一个迭代data.frame对象的函数。 它只是添加一个新的列data.frame和积累的东西。 （操作简单）。 data.frame大约有850K行。 我的电脑还在工作（现在大约10小时），我不知道运行时间。
dayloop2 <- function(temp){ for (i in 1:nrow(temp)){ temp[i,10] <- i if (i > 1) { if ((temp[i,6] == temp[i-1,6]) & (temp[i,3] == temp[i-1,3])) { temp[i,10] <- temp[i,9] + temp[i-1,10] } else { temp[i,10] <- temp[i,9] } } else { temp[i,10] <- temp[i,9] } } names(temp)[names(temp) == "V10"] <- "Kumm." return(temp) }
任何想法如何加快这一行动？
最大的问题和无效的根源是索引data.frame，我的意思是所有这些行你使用temp[,] 。
尽量避免这种情况。 我把你的function，改变索引，在这里version_A
dayloop2_A <- function(temp){ res <- numeric(nrow(temp)) for (i in 1:nrow(temp)){ res[i] <- i if (i > 1) { if ((temp[i,6] == temp[i-1,6]) & (temp[i,3] == temp[i-1,3])) { res[i] <- temp[i,9] + res[i-1] } else { res[i] <- temp[i,9] } } else { res[i] <- temp[i,9] } } temp$`Kumm.` <- res return(temp) }
正如你所看到的，我创build了收集结果的vectorres 。 最后我把它添加到data.frame ，我不需要data.frame名称。 那么它有多好？
我运行每个函数的data.frame从1,000到10,000 1,000，用system.time测量时间
X <- as.data.frame(matrix(sample(1:10, n*9, TRUE), n, 9)) system.time(dayloop2(X))
结果是
你可以看到你的版本取决于nrow(X)指数。 修改后的版本具有线性关系，而简单的lm模型预测对于85万行计算需要6分10秒。
vector化的力量
正如Shane和Calimo所说的，向量化是提高性能的关键。 从你的代码中你可以移出循环：
- 空调
- 结果的初始化（
temp[i,9]）
这导致这个代码
dayloop2_B <- function(temp){ cond <- c(FALSE, (temp[-nrow(temp),6] == temp[-1,6]) & (temp[-nrow(temp),3] == temp[-1,3])) res <- temp[,9] for (i in 1:nrow(temp)) { if (cond[i]) res[i] <- temp[i,9] + res[i-1] } temp$`Kumm.` <- res return(temp) }
比较这个function的结果，这次是10000到nrow 。
调整调谐
另一个调整是将循环索引temp[i,9]为res[i] （这在第i次循环迭代中是完全相同的）。 索引一个向量和索引一个数据data.frame又是不同的。
第二件事：当你看循环时，你可以看到没有必要循环所有的i ，但只适合那些符合条件的。
所以我们走了
dayloop2_D <- function(temp){ cond <- c(FALSE, (temp[-nrow(temp),6] == temp[-1,6]) & (temp[-nrow(temp),3] == temp[-1,3])) res <- temp[,9] for (i in (1:nrow(temp))[cond]) { res[i] <- res[i] + res[i-1] } temp$`Kumm.` <- res return(temp) }
你获得的性能高度取决于数据结构。 准确地说 – 在条件中TRUE值的百分比。 对于我的模拟数据，需要一秒以下的时间计算850,000行。
我想你可以走得更远，至less有两件事情可以做到：
- 写一个
C代码做条件cumsum -
如果你知道在你的数据最大序列不是很大，那么你可以改变循环向量化，而类似的东西
while (any(cond)) { indx <- c(FALSE, cond[-1] & !cond[-n]) res[indx] <- res[indx] + res[which(indx)-1] cond[indx] <- FALSE }
用于模拟和数字的代码可以在GitHub上find 。
加速R代码的一般策略
首先，找出慢速的部分。 没有必要优化运行缓慢的代码。 对于less量的代码，只需简单的思考就可以工作。 如果失败了，RProf和类似的分析工具可能会有所帮助。
一旦你找出了瓶颈，想想更有效的algorithm来做你想做的事情。 如果可能的话，计算应该只运行一次，所以：
- 存储结果并访问它们而不是反复重新计算
- 循环中取决于非循环的计算
- 避免不必要的计算（例如，不要使用固定search的正则expression式 ）
使用更高效的function可以产生中等或大的速度增益。 例如， paste0产生一个小的效率增益，但是.colSums()和它的亲戚产生更明显的增益。 mean是特别慢 。
那么你可以避免一些特别常见的麻烦 ：
-
cbind会很快让你放慢脚步。 - 初始化您的数据结构，然后填充它们， 而不是每次都展开它们 。
- 即使预先分配，您也可以切换到逐个引用的方式，而不是按值传递的方式，但这可能不值得一提。
- 看看R地狱更多的陷阱，以避免。
尝试更好的vector化 ，这往往可以但不总是帮助。 在这方面，像ifelse ， diff等内在向量化的命令将比apply的命令系列（在很好的写入循环中几乎不提供速度提升）提供更多的改进。
您也可以尝试向Rfunction提供更多信息 。 例如，使用vapply而不是sapply ，并且在基于文本的数据中读取时指定colClasses 。 速度的提高将取决于你消除了多less猜测。
接下来，考虑优化软件包 ： data.table软件包可以在数据操作和读取大量数据（ fread ）的情况下产生巨大的速度收益。
接下来，尝试通过更高效的方式调用R来获得更快的速度：
- 编译你的R脚本。 或者使用
Ra和jit包进行即时编译（Dirk在本演示中有一个例子）。 - 确保你正在使用优化的BLAS。 这些提供了全面的速度收益。 老实说，R不能自动使用最有效的安装库。 希望Revolution R将他们在这里所做的工作贡献给整个社区。
- Radford Neal做了一系列优化，其中一些被优化，R Core被采用，还有很多被分成了pqR 。
最后，如果以上所有内容仍然不能满足您的需求，那么您可能需要为较慢的代码片段使用更快的语言 。 Rcpp和inline在这里的结合使得用C ++代码只replacealgorithm的最慢部分特别容易。 例如，这是我第一次尝试这样做 ，它甚至吹捧了高度优化的R解决scheme。
如果你仍然留下麻烦，你只需要更多的计算能力。 看看并行化 （ http://cran.r-project.org/web/views/HighPerformanceComputing.html ），甚至是基于GPU的解决scheme（ gpu-tools ）。
链接到其他指导
如果你正在使用for循环，你很可能编码R，就像是C或Java或其他。 正确vector化的R代码速度非常快。
以这两个简单的代码位为例，按顺序生成一个10000个整数的列表：
第一个代码示例是如何使用传统的编码范例编码循环。 完成需要28秒
system.time({ a <- NULL for(i in 1:1e5)a[i] <- i }) user system elapsed 28.36 0.07 28.61
通过预分配内存的简单操作，您可以获得几乎100倍的提升：
system.time({ a <- rep(1, 1e5) for(i in 1:1e5)a[i] <- i }) user system elapsed 0.30 0.00 0.29
但是使用使用冒号运算符的基本R向量操作:这个操作实际上是瞬时的：
system.time(a <- 1:1e5) user system elapsed 0 0 0
通过使用索引或嵌套的ifelse()语句跳过循环可以使速度更快。
idx <- 1:nrow(temp) temp[,10] <- idx idx1 <- c(FALSE, (temp[-nrow(temp),6] == temp[-1,6]) & (temp[-nrow(temp),3] == temp[-1,3])) temp[idx1,10] <- temp[idx1,9] + temp[which(idx1)-1,10] temp[!idx1,10] <- temp[!idx1,9] temp[1,10] <- temp[1,9] names(temp)[names(temp) == "V10"] <- "Kumm."
正如Ari在回答结束时提到的那样， Rcpp和inline软件包让事情变得非常容易。 作为一个例子，试试这个inline代码（警告：未testing）：
body <- 'Rcpp::NumericMatrix nm(temp); int nrtemp = Rccp::as<int>(nrt); for (int i = 0; i < nrtemp; ++i) { temp(i, 9) = i if (i > 1) { if ((temp(i, 5) == temp(i - 1, 5) && temp(i, 2) == temp(i - 1, 2) { temp(i, 9) = temp(i, 8) + temp(i - 1, 9) } else { temp(i, 9) = temp(i, 8) } } else { temp(i, 9) = temp(i, 8) } return Rcpp::wrap(nm); ' settings <- getPlugin("Rcpp") # settings$env$PKG_CXXFLAGS <- paste("-I", getwd(), sep="") if you want to inc files in wd dayloop <- cxxfunction(signature(nrt="numeric", temp="numeric"), body-body, plugin="Rcpp", settings=settings, cppargs="-I/usr/include") dayloop2 <- function(temp) { # extract a numeric matrix from temp, put it in tmp nc <- ncol(temp) nm <- dayloop(nc, temp) names(temp)[names(temp) == "V10"] <- "Kumm." return(temp) }
#include东西有一个类似的过程，你只要传递一个参数即可
inc <- '#include <header.h>
cxxfunction，如include=inc 。 这真的很酷，它为你做了所有的链接和编译，所以原型是非常快的。
免责声明：我不完全确定tmp的类应该是数字，而不是数字matrix或其他。 但我很确定。
编辑：如果在这之后你仍然需要更多的速度， OpenMP是一个适用于C++的并行化工具。 我没有尝试从inline使用它，但它应该工作。 这个想法是，在n核的情况下，循环迭代k由k % n来执行。 在Matloff的“艺术的R编程”一书中可以find一个合适的介绍，在这里可以在第16章中find。
我不喜欢重写代码…当然ifelse和lapply是更好的select，但有时很难做到这一点。
我经常使用data.frames作为使用列表，如df$var[i]
这是一个制造的例子：
nrow=function(x){ ##required as I use nrow at times. if(class(x)=='list') { length(x[[names(x)[1]]]) }else{ base::nrow(x) } } system.time({ d=data.frame(seq=1:10000,r=rnorm(10000)) d$foo=d$r d$seq=1:5 mark=NA for(i in 1:nrow(d)){ if(d$seq[i]==1) mark=d$r[i] d$foo[i]=mark } }) system.time({ d=data.frame(seq=1:10000,r=rnorm(10000)) d$foo=d$r d$seq=1:5 d=as.list(d) #become a list mark=NA for(i in 1:nrow(d)){ if(d$seq[i]==1) mark=d$r[i] d$foo[i]=mark } d=as.data.frame(d) #revert back to data.frame })
data.frame版本：
user system elapsed 0.53 0.00 0.53
列表版本：
user system elapsed 0.04 0.00 0.03
使用vector列表的速度是data.frame的17倍。
有关为什么内部data.frames在这方面如此缓慢的任何意见？ 人们会认为他们像列表一样运作
为了更快的代码，这个class(d)='list'而不是d=as.list(d)和class(d)='data.frame'
system.time({ d=data.frame(seq=1:10000,r=rnorm(10000)) d$foo=d$r d$seq=1:5 class(d)='list' mark=NA for(i in 1:nrow(d)){ if(d$seq[i]==1) mark=d$r[i] d$foo[i]=mark } class(d)='data.frame' }) head(d)
在R中，通常可以使用apply系列函数加速循环处理（在你的情况下，它可能会被replicate ）。 看看提供进度条的plyr包。
另一种select是完全避免循环，并用vector化算术replace它们。 我不确定你到底在做什么，但是你可以把你的函数同时应用到所有的行：
temp[1:nrow(temp), 10] <- temp[1:nrow(temp), 9] + temp[0:(nrow(temp)-1), 10]
这将会快得多，然后你可以用你的条件过滤行：
cond.i <- (temp[i, 6] == temp[i-1, 6]) & (temp[i, 3] == temp[i-1, 3]) temp[cond.i, 10] <- temp[cond.i, 9]
向量化的算术需要更多的时间和思考这个问题，但是有时你可以在执行时间内节省几个数量级。
用data.table处理是一个可行的select：
n <- 1000000 df <- as.data.frame(matrix(sample(1:10, n*9, TRUE), n, 9)) colnames(df) <- paste("col", 1:9, sep = "") library(data.table) dayloop2.dt <- function(df) { dt <- data.table(df) dt[, Kumm. := { res <- .I; ifelse (res > 1, ifelse ((col6 == shift(col6, fill = 0)) & (col3 == shift(col3, fill = 0)) , res <- col9 + shift(res) , # else res <- col9 ) , # else res <- col9 ) } ,] res <- data.frame(dt) return (res) } res <- dayloop2.dt(df) m <- microbenchmark(dayloop2.dt(df), times = 100) #Unit: milliseconds # expr min lq mean median uq max neval #dayloop2.dt(df) 436.4467 441.02076 578.7126 503.9874 575.9534 966.1042 10
如果您忽略条件过滤的可能收益，则速度非常快。 很显然，如果你可以对数据的子集进行计算，那么它是有帮助的。
