R - 需要帮助加速 for 循环

Question

R - 需要帮助加速 for 循环

rfor-loopdataframevectorizationdata.table

5

我有两个数据框; 其中一个包含48行，看起来像这样:

名称 = Z31

     Est.Date   Site    Cultivar   Planting
1   24/07/2011  Birchip Axe           1
2   08/08/2011  Birchip Bolac         1
3   24/07/2011  Birchip Derrimut      1
4   12/08/2011  Birchip Eaglehawk     1
5   29/07/2011  Birchip Gregory       1
6   29/07/2011  Birchip Lincoln       1
7   23/07/2011  Birchip Mace          1
8   29/07/2011  Birchip Scout         1
9   17/09/2011  Birchip Axe           2
10  19/09/2011  Birchip Bolac         2

另一个文件包含模拟器的输出，行数超过23000行。它的格式如下所示：

名称 = 预测

    Date        maxt    mint    Cultivar    Site    Planting    tt  cum_tt
1   5/05/2011   18       6.5    Axe        Birchip  1        12.25  0
2   6/05/2011   17.5     2.5    Axe        Birchip  1        10     0
3   7/05/2011   18       2.5    Axe        Birchip  1        10.25  0
4   8/05/2011   19.5       2    Axe        Birchip  1        10.75  0
5   9/05/2011   17       4.5    Axe        Birchip  1        10.75  0
6   10/05/2011  15.5    -0.5    Axe        Birchip  1        7.5    0
7   11/05/2011  14       5.5    Axe        Birchip  1        9.75   0
8   12/05/2011  19         8    Axe        Birchip  1        13.5   0
9   13/05/2011  18.5     7.5    Axe        Birchip  1        13     0
10  14/05/2011  16       3.5    Axe        Birchip  1        9.75   0

我想要做的是让cum_tt列开始将当前行的tt列加到前一行的cum_tt中（累积相加），仅当pred DF中的日期等于或大于Z31 est.Date时。我已经编写了以下for循环：

for (i in 1:nrow(Z31)){
  for (j in 1:nrow(pred)){
    if (Z31[i,]$Site == pred[j,]$Site & Z31[i,]$Cultivar == pred[j,]$Cultivar & Z31[i,]$Planting == pred[j,]$Planting &
        pred[j,]$Date >= Z31[i,]$Est.Date)
    {
      pred[j,]$cum_tt <- pred[j,]$tt + pred[j-1,]$cum_tt
    }
  }
}

这段代码可以运行，但速度太慢，需要花费大约一个小时才能完成整个操作。我知道循环不是R语言的强项，所以能否有人帮助我将此操作向量化?

提前致谢。

更新

以下是dput（Z31）的输出：

结构(list(Est.Date = 结构(c(15179, 15194, 15179, 15198, 15184, 15184, 15178, 15184, 15234, 15236, 15230, 15238, 15229, 15236, 15229, 15231, 15155, 15170, 15160, 15168, 15165, 15159, 15170, 15170, 15191, 15205, 15198, 15203, 15202, 15195, 15203, 15206, 15193, 15193, 15195, 15200, 15193, 15205, 15200, 15205, 15226, 15245, 15231, 15259, 15241, 15241, 15241, 15241), class = "Date"), Site = 结构(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), .Label = 结构(c("Birchip", "Gatton", "Tarlee"), class = "factor"), class = "factor"), Cultivar = 结构(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L), .Label = 结构(c("Axe", "Bolac", "Derrimut", "Eaglehawk", "Gregory", "Lincoln", "Mace", "Scout"), class = "factor"), class = "factor"), Planting = c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L)), .Names = 结构(c("Est.Date", "Site", "Cultivar", "Planting"), row.names = c(NA, -48L), class = "data.frame")

这是pred。注意这里有一些额外的列。为了方便阅读，我只包含了相关的列。

structure(list(Date = structure(c(15099, 15100, 15101, 15102, 
15103, 15104, 15105, 15106, 15107, 15108, 15109, 15110, 15111, 
15112, 15113, 15114, 15115, 15116, 15117, 15118), class = "Date"), 
    flowering_das = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 
    0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L), Zadok = c(9, 9, 
    9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 11, 11.032, 11.085, 
    11.157), stagename = structure(c(8L, 8L, 8L, 8L, 8L, 8L, 
    8L, 8L, 9L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 3L, 3L, 3L), .Label = c("emergence", 
    "end_grain_fill", "end_of_juvenil", "floral_initiat", "flowering", 
    "germination", "maturity", "out", "sowing", "start_grain_fi"
    ), class = "factor"), node_no = c(0, 0, 0, 0, 0, 0, 0, 0, 
    0, 0, 0, 0, 0, 0, 0, 0, 2, 2.032, 2.085, 2.157), maxt = c(18, 
    17.5, 18, 19.5, 17, 15.5, 14, 19, 18.5, 16, 16, 15, 16.5, 
    16.5, 20.5, 23, 25.5, 16.5, 16.5, 15), mint = c(6.5, 2.5, 
    2.5, 2, 4.5, -0.5, 5.5, 8, 7.5, 3.5, 6, 1, 5.5, 2, 7, 7, 
    9, 13.5, 11.5, 8.5), Cultivar = c("Axe", "Axe", "Axe", "Axe", 
    "Axe", "Axe", "Axe", "Axe", "Axe", "Axe", "Axe", "Axe", "Axe", 
    "Axe", "Axe", "Axe", "Axe", "Axe", "Axe", "Axe"), Site = c("Birchip", 
    "Birchip", "Birchip", "Birchip", "Birchip", "Birchip", "Birchip", 
    "Birchip", "Birchip", "Birchip", "Birchip", "Birchip", "Birchip", 
    "Birchip", "Birchip", "Birchip", "Birchip", "Birchip", "Birchip", 
    "Birchip"), Planting = c("1", "1", "1", "1", "1", "1", "1", 
    "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", 
    "1"), `NA` = c("Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", 
    "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", 
    "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", 
    "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", 
    "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", 
    "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", 
    "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out", "Birchip TOS1 Axe.out"
    ), tt = c(12.25, 10, 10.25, 10.75, 10.75, 7.5, 9.75, 13.5, 
    13, 9.75, 11, 8, 11, 9.25, 13.75, 15, 17.25, 15, 14, 11.75
    ), cum_tt = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
    0, 0, 0, 0, 0)), .Names = c("Date", "flowering_das", "Zadok", 
"stagename", "node_no", "maxt", "mint", "Cultivar", "Site", "Planting", 
NA, "tt", "cum_tt"), row.names = c(NA, 20L), class = "data.frame")

更新

感谢大家的帮助。我对向量方式的操作还不太熟悉，所以没能及时实现一些更复杂的解决方案。下面是按Subs建议的方法计时的结果。现在速度已经足够快，可以满足我的需求了。以下数字表示Z在P上进行一次迭代所需的时间（单位为秒）。

我的方式：59.77

Subs的方式：14.62

使用数值日期的Subs方式：11.12

- Justin

请问您能否更新您的问题，附上dput(Z31)和dput(pred)的结果？在阅读文本后，我无法按原样运行您的代码... - Chase

1

顺便提一下，使用data.table()和roll = TRUE看起来会更快。 - Chase

关于dput输出，与上述相同。看了一下data.table...午饭后可以探索一下 :) - Justin

dput(head(pred,20)) 应该可以解决问题 :) - Chase

当然是 bonk 了... 给你 :) 现在去吃饭了。 - Justin

如果你想比较日期，为什么不将它们转换为向量而不是“日期类”？这样很容易进行排序和比较。 - Subs

3个回答

1

一个快速的解决方案是在循环外面定义一个向量，如下所示：

 temp_cumtt=c(rep(0,nrow(pred)))

然后使用这个：

if (Z31[i,2] == pred[j,5] & Z31[i,3] == pred[j,4] & Z31[i,4] == pred[j,6] & pred[j,1] >= Z31[i,1]){
   temp_cumtt[j]=pred[j,7] + pred[j-1,8]
}

不要直接更新数据框列。

当你退出循环后，可以更新该列：

 pred$cum_tt = temp_cumtt

另一件事是，在使用索引为1的j时，要小心使用j-1。在您的示例中，它不会导致条件问题。

编辑：

现在看着您的数据格式，我有以下建议。

1）不要转换为Date类，而是将其保留为值向量。

2）根据Date向量对Z31 data_frame进行排序：Z31=Z31[with(Z31, order(-Date)), ]（请注意，按降序排列，因为您想比较pred [，Date index] >= Z31 [，Date index]）

3）将第一个循环用作pred。首先取pred的日期->pred[i,1]，然后尝试二进制排序并找到满足Z31的哪个索引，并从该索引开始向下移动列表。如果条件满足，则检查其余条件并像以前一样填充temp_cumtt[i]。

这应该是非常快的（因为二进制排序仅在48行Z31上进行，您可以将运行时间与其他解决方案进行比较）。

- Subs

感谢订阅，这确实提高了执行速度，我从每次外部迭代需要60秒的速度降到了约15秒。虽然循环非常慢，但我希望能够完全不使用循环来完成这个任务。尽管如此，这是一个巨大的改进。我不担心 j-1 索引，因为在我使用的数据集中，第一行永远不会成立。 - Justin

循环并不低效，它取决于你在其中做什么。另一种方法是根据你的日期对Z31数据框进行排序，由于pred中的日期已经排序，因此你可以使用算法从一个索引开始比较，而不是遍历所有的Z行* P行。 - Subs

谢谢你的帮助，subs。我需要继续我的工作，所以我无法实现二进制排序技术。不过你已经帮助我将其速度提高到了可以接受的水平。以上问题中有一些时间限制。 - Justin

0

让我们利用 data.table，这应该大大加快速度。

Z31 <- data.table(Z31,key="Site,Cultivar,Planting")
pred <- data.table(pred)

## First, let's create an extra column in `pred` to see the corresponding date from `Z31` 
## Note 1: The JT is necessary since both sets have the same column names
## Note 2: I needed to use as.integer on Planting to make it work

pred[,Z31Est.Date:={JT=J(Site,Cultivar,as.integer(Planting)); Z31[JT,Est.Date][[4]]}]

## Now we can see for each row whether the date in `pred` is higher than or equal to that from `Z31`.

pred[,DateTrue:=Date>=Z31Est.Date]

## Finally, we only have to add up `pred[i,tt]` and `pred[i-1,cum_tt]` for each row where `DateTrue` equals `TRUE`.

for (i in 1:nrow(pred)) set(pred,i,13L,if(pred[i,DateTrue]) pred[i-1,cum_tt]+pred[i,tt] else(0))

- Dirk

一个仍然包含 for (i in 1:nrow(pred)) 的“解决方案”并不是真正有用的——我们已经有了一个完全向量化的解决方案。 - Tim P

问题是要加速循环，对吧？我认为这样可以加速它，但我无法测试它是否真的有效。但确实应该避免使用循环。 - Dirk

我理解这个问题是...“我知道循环不是R的强项，所以有谁能帮助我向量化这个操作？” :) - Tim P

@Dirk 我能看出你在这里的尝试。虽然我怀疑可能存在一个更短、更简单的 data.table 解决方案。甚至可能只需要一两行代码，比 Tim 使用 paste 进行向量化的方法更快（通常是无法扩展到 1e6、1e7 或 1e9 行的部分）。我会回来处理这个问题的。我们可以使用 1e7 的示例来演示时间。 - Matt Dowle

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Tim P · Accepted Answer

当然，我们可以在几秒钟内完成这个...这是我在这里的第一个答案，请温柔一点！

## first make sure we have dates in a suitable format for comparison
## by using strptime, creating the columns estdate_tidy and date_tidy
## in Z31 and pred respectively

Z31$estdate_tidy = strptime(as.character(Z31$Est.Date), "%d/%m/%Y")
pred$date_tidy = strptime(as.character(pred$Date), "%d/%m/%Y")

## now map the estdate_tidy column over to pred using the match command -
## Z31_m and pred_m are dummy variables that hopefully make this clear

Z31_m = paste(Z31$Site, Z31$Cultivar, Z31$Planting)
pred_m = paste(pred$Site, pred$Cultivar, pred$Planting)
pred$estdate_tidy = Z31$estdate_tidy[match(pred_m, Z31_m)]

## then define a ttfilter column that copies tt, except for being 0 when
## estdate_tidy is after date_tidy (think this is what you described)

pred$ttfilter = ifelse(pred$date_tidy >= pred$estdate_tidy, pred$tt, 0)

## finally use cumsum function to sum this new column up (looks like you
## wanted the answer in cum_tt so I've put it there)

pred$cum_tt = cumsum(pred$ttfilter)

希望这可以帮到你 :)

更新（6月7日）：

用于处理新规范的矢量化代码 - 即应分别针对每组条件（站点/品种/种植）进行累加 - 如下所示：

Z31$Lookup=with(Z31,paste(Site,Cultivar,Planting,sep="~"))
Z31$LookupNum=match(Z31$Lookup,unique(Z31$Lookup))
pred$Lookup=with(pred,paste(Site,Cultivar,Planting,sep="~"))
pred$LookupNum=match(pred$Lookup,unique(pred$Lookup))

pred$Est.Date = Z31$Est.Date[match(pred$Lookup,Z31$Lookup)]
pred$ttValid = (pred$Date>=pred$Est.Date)
pred$ttFiltered = ifelse(pred$ttValid, pred$tt, 0)

### now fill in cumsum of ttFiltered separately for each LookupNum
pred$cum_tt_Z31 = as.vector(unlist(tapply(pred$ttFiltered,
                                          pred$LookupNum,cumsum)))

我的电脑运行时间为0.16秒，最终的pred$cum_tt_Z31列与非矢量化代码的答案完全匹配 :)

为了完整起见，值得注意的是，上述最终复杂的tapply行可以用以下更简单的方法替换，只需在48种可能情况下进行简短的循环：

pred$cum_tt_Z31 = rep(NA, nrow(pred))
for (lookup in unique(pred$Lookup)) {
    subs = which(pred$Lookup==lookup)
    pred$cum_tt_Z31[subs] = cumsum(pred$ttFiltered[subs])
}

由于这里的循环非常小且循环内的工作已经向量化，因此运行时间仅略微增加到0.25秒左右。

我想我们解决了它！ :)

有关向量化的一些快速观察（6月8日）：

将过程步骤向量化后，总运行时间从近1小时降至0.16秒。即使考虑不同机器速度，这也是至少比通过进行微小调整但仍保留循环结构获得的2-5倍速度提升要高出10,000倍，这远远超过了前者。

第一个关键观察：在解决方案中，每行都会创建一个与Z31或pred列长度相同的全新向量，而无需循环。为了整洁，我通常发现将这些新向量创建为新的数据框列很有用，但显然这并非必须的。

第二个观察：使用“粘贴并匹配”策略将所需的Est.Date列从Z31正确地传输到pred。对于这种任务，还有其他方法（例如使用merge），但我选择这条路线，因为它完全可靠，并保证保留pred中的顺序（这是至关重要的）。基本上，粘贴操作只是让您一次匹配多个字段，因为如果粘贴的字符串匹配，则它们所有的组成部分都匹配。我使用~作为分隔符（前提是我知道该符号不会出现在任何字段中），以避免粘贴操作引起任何歧义。如果您使用空格分隔符，则将("A B"，"C"，"D")粘合在一起会给出与("A"，"B C"，"D")粘合在一起相同的结果-我们要避免任何麻烦！

第三个观察：很容易将逻辑运算向量化，例如查看一个向量是否超过另一个向量（请参见pred$ttValid），或根据向量的值选择一个值或另一个值（请参见pred$ttFiltered）。在当前情况下，这些可以合并为一行，但我将其拆分为更详细的演示。

第四个观察：最后一行创建了pred$cum_tt_Z31，基本上只是通过tapply沿着对应于每个单独的pred$LookupNum值的行执行cumsum操作，tapply允许您在不同的行组上应用相同的函数（在这里，我们按pred$LookupNum进行分组）。pred$LookupNum的定义在这里非常有帮助-它是一个具有一块1、一块2等的数字索引。这意味着tapply产生的cumsum向量的最终列表可以简单地取消列表并放入向量中，自动处于正确的顺序。如果您通过未按此方式排序的组进行tapply和拆分，则通常需要几行额外的代码才能正确地将事物匹配回来（尽管这并不难）。

最终观察：如果最后的tapply让人感到害怕，值得强调的是，如果循环中的工作被很好地向量化，则在少量案例（例如48个）上快速循环并不总是会导致灾难。更新部分末尾的“替代方法”显示，可以通过预先准备一个答案列（最初全部为NA），然后逐个处理48个子集并使用适当的cumsum填充每个块来实现组累加步骤。正如文本中所述，这一单步大约只有聪明的tapply方法的一半速度，并且如果需要更多的子集，这将变得相当慢。

如果有任何关于这种任务的后续问题，请不要犹豫给我发消息。