在R (Rcpp)中提高4嵌套for循环的速度？
让我们假设我有一个包含1000个条目/行的数据帧。每一行都有一个ID，第二列包含一些数据，第三列也包含一些数据。
所以数据框看起来像这样：
ID yesNo Id_specific_data
1 1 4
2 0 8
3 0 43
4 1 11
5 0 9..。诸若此类。
我现在需要执行以下操作：
n = 4
ID_range <- c(1:n)
ID_spec_data <- floor(runif(n, min=10, max=100))
yesNo_data <- sample(c(0,1), replace=TRUE, size=n)
df <- data.frame("ID" = ID_range, "yesNo" = yesNo_data, "ID_specific_data" = ID_spec_data)
m <- 1
for (i in seq(1, 100, 1)) {
for (j in seq(0.1, 1, 0.1)) {
log_like_list <- c()
for (k in seq(0.1, 1, 0.1)) {
total_ID_list <- c()
for (l in seq(1, length(df$ID))) {
x = (df$ID_specific_data[[l]]*k - j) / (i*j)
calc = pnorm(x, mean = 0, sd = 1, lower.tail = TRUE, log.p = FALSE)
total_ID_list[[l]] = calc
}
# log likelihood function
final_calc = sum(df$yesNo*log(total_ID_list)+(1-df$yesNo)*log(1 - total_ID_list))
log_like_list[[m]] = final_calc
m <- m + 1
}
}
}因此，基本上最终结果(log_like_list)应该是一个带有1500*200*100值的列表/向量。但为了做到这一点，需要对数据帧中ID的数量进行相同数量的计算(在我的情况下大约是500-1000 )。总而言之--大量的计算。
我知道就速度而言，for循环可能是你能做的最糟糕的事情，但我甚至不确定在计算如此之多的情况下，使用apply是否会让它变得超级快？我读过关于Rcpp的文章，它原则上可以减少计算时间，是所有选项中最重要的。但它需要我所看到的C++知识(这是我真正缺乏的)，我甚至不确定它是否适用于我的问题？
那么，是否可以通过任何R技巧显著减少计算时间，或者我只需等待它结束？
回答 2
Stack Overflow用户
发布于 2018-08-27 05:28:23
我认为您当前的编辑仍然是错误的，您可能不应该在任何循环中重新定义log_like_list。这里有一个替代方案，首先使用expand.grid分配所有参数组合，这在内存方面有点浪费，但我认为它是可管理的：
n <- 4L
df <- data.frame(
ID = 1L:n,
yesNo = sample(c(0,1), replace=TRUE, size=n),
ID_specific_data = floor(runif(n, min=10, max=100))
)
params <- expand.grid(
i = seq(1, 100, 1),
j = seq(0.1, 1, 0.1),
k = seq(0.1, 1, 0.1)
)
log_like <- sapply(1L:nrow(params), function(row_id) {
i <- params$i[row_id]
j <- params$j[row_id]
k <- params$k[row_id]
calc <- sapply(df$ID_specific_data, function(idsd) {
x <- (idsd * k - j) / (i * j)
pnorm(x, mean = 0, sd = 1, lower.tail = TRUE, log.p = FALSE)
})
sum(df$yesNo * log(calc) + (1 - df$yesNo) * log(1 - calc))
})但是，对于您的最终用例，这可能仍然太慢了……您可以尝试使用并行化，如果您有多个内核，则可能会有可接受的时间：
library(doParallel)
library(itertools)
# do NOT run these lines several times without calling stopCluster() on the created workers
workers <- makeCluster(detectCores())
registerDoParallel(workers)
n <- 1000L
df <- data.frame(
ID = 1L:n,
yesNo = sample(c(0,1), replace=TRUE, size=n),
ID_specific_data = floor(runif(n, min=10, max=100))
)
params <- expand.grid(
i = seq(1, 150, 0.1),
j = seq(0.1, 2, 0.01),
k = seq(0.1, 1, 0.01)
)
params_chunk <- isplitRows(params, chunks = getDoParWorkers())
log_like_par <- foreach(param = params_chunk, .combine = c, .multicombine = TRUE) %dopar% {
# return from foreach body here
sapply(1L:nrow(param), function(row_id) {
i <- param$i[row_id]
j <- param$j[row_id]
k <- param$k[row_id]
calc <- sapply(df$ID_specific_data, function(idsd) {
x <- (idsd * k - j) / (i * j)
pnorm(x, mean = 0, sd = 1, lower.tail = TRUE, log.p = FALSE)
})
# return from sapply body here
sum(df$yesNo * log(calc) + (1 - df$yesNo) * log(1 - calc))
})
}
stopCluster(workers); registerDoSEQ()我试着在我的系统(4个内核)中运行它，但几分钟后就停止了。如果你等它结束，让我知道它花了多长时间。
Stack Overflow用户
发布于 2018-08-26 22:34:05
这不会是一个100%的答案，你可以复制和粘贴，但我认为它将帮助您获得部分方法。主要是，你需要思考为什么你要花时间做循环，而你实际上是在处理本质上不变的值。
例如
i <- seq(1, 100, 1)
j <- seq(0.1, 1, 0.1)
ioxj <- i %o% j
df_ij <- data.frame("i" = i, "j" = j, "ioxj" = ioxj)
df_ij$ixj <- df_ij$i * df_ij$j将得到i和j的所有组合以及它们的乘积，没有理由使用循环来获得基本的数学结果。您可以在某个时刻使用循环遍历这些列，这可能是有意义的，因为i和j的值可能会改变。您也可以使用k进行类似的操作。
此外，在循环中遍历数据帧中的每一行，也没有理由像这样做x = (df$ID_specific_data[[l]]*k - j) / (i*j)，这会失去矢量化的整个概念，你想要以此来结束。x = (df$ID_specific_data*k - j) / (i*j)
您需要摆弄代码以获得您想要的方式，但花时间这样做是值得的。偶尔的循环可能是正确的，但我认为你最终可能会做一些简单得多的事情。
https://stackoverflow.com/questions/52022473
复制相似问题