如何在R中高效地扩展向量?

3
假设您有一个函数,它以数字作为输入并输出向量。然而,输出向量的大小取决于输入,并且在函数之前无法计算。
例如,考虑 3N+1著名算法。该算法的简单实现,返回到1的整个路径可能如下所示:
compute <- function(x) {
  if (x %% 2 == 0)
    return(x / 2)
  return(3*x + 1)
}

algo <- function(x) {
  if (x == 1)
    return(1)

  output <- x
  while(x != 1) {
    x <- compute(x)
    output <- c(output, x)
  }

  return(output)
}
algo函数根据函数的要求返回输入X到1的整个路径。正如您所看到的,output变量使用c()(合并)函数动态增长。

是否有其他替代方法?增加列表速度更快吗?我应该采用一些经典的动态向量逻辑,比如初始化一个空的N大小的向量,并在每次满时将其加倍?

编辑:请不要试图优化我的辅助函数的结构。我明白,但这不是重点!我只关心c()函数及其替代方案。


转向分配多少内存,对于初始值高达1,000,000,使用n <- sapply(1:1e6, function(x) length(algo(x))),输出向量的最大长度为525。 - dcarlson
请注意,这只是一个例子,我并不特别想解决这个问题。这只是一个小故事! - eduardokapp
1
请查看 https://privefl.github.io/blog/why-loops-are-slow-in-r/。 - F. Privé
@F.Privé 有趣的文章,但我想要添加额外的警告,即 JIT 编译在 R 中并不是通常所指的 JIT。也就是说,它没有编译成本地代码。 - user13963867
3个回答

1

更新

根据您的编辑,也许您可以检查以下解决方案

algo_TIC2 <- function(x) {
  res <- x
  repeat {
    u <- tail(res, 1)
    if (u != 1) {
      res[length(res) + 1] <- if (u %% 2) 3 * u + 1 else u / 2
    } else {
      return(res)
    }
  }
}

您可以像下面这样使用递归:
compute <- function(x) if (x %% 2) 3*x + 1 else x / 2
algo_TIC1 <- function(x) {
  if (x == 1) {
    return(1)
  }
  c(x, algo_TIC1(compute(x)))
}

并且你会看到

> algo_TIC1(3000)
 [1] 3000 1500  750  375 1126  563 1690  845 2536 1268  634  317  952  476  238
[16]  119  358  179  538  269  808  404  202  101  304  152   76   38   19   58
[31]   29   88   44   22   11   34   17   52   26   13   40   20   10    5   16
[46]    8    4    2    1

如果您不想使用任何辅助函数,例如compute,您可以尝试。
algo_TIC1 <- function(x) {
  if (x == 1) {
    return(1)
  }
  c(x, algo_TIC1(if (x %% 2) 3*x + 1 else x / 2))
}

这很优雅,但它是否解决了他对于向量大小的“效率担忧”? - Dirk Eddelbuettel
@DirkEddelbuettel 看起来好像不是 :( - ThomasIsCoding
在这里使用ifelse()没有意义。只需使用if (x %% 2) 3*x + 1 else x / 2。对于计算的这一部分,它大约快了4倍。 - user2554330
@user2554330 非常感谢您的建议。确实更快了! - ThomasIsCoding
谢谢您的建议,但我不认为这解决了我的问题。实际上,就动态增长向量而言,它看起来几乎一样。 - eduardokapp

1

所以,困扰你的是重新分配,你是正确的。让我们来看看。

library(microbenchmark)

microbenchmark({
  a <- c()
  for (i in seq(1e4)) {
    a <- c(a, i)
  }
})

microbenchmark({
  a <- numeric(1e4)
  for (i in seq(1e4)) {
    a[[i]] <- i
  }
})

microbenchmark({
  a <- numeric(1)
  k <- 1
  for (i in seq(1e4)) {
    if (i > k) {
      a <- c(a, numeric(k))
      k <- k + k
    }
    a[[i]] <- i
  }
  a <- head(a, 1e4)
})

而且时间如下:

Append
     min       lq      mean   median       uq      max neval
 78.0162 78.67925  83.36224 79.54515 81.79535 166.6988   100

Preallocate
     min       lq     mean    median       uq      max neval
1.484901 1.516051 1.567897    1.5552 1.569451 1.895601   100

Amortize
     min       lq     mean    median       uq      max neval
3.316501 3.377201  3.62415  3.484351 3.585701  11.7596   100

不要向向量中追加太多元素。如果可能的话,预先分配内存,否则就进行摊销分配。

即使你事先不知道实际大小,也可能有一个上限。然后你仍然可以预先分配内存,并在最后截断。即使是一个合理的估计也是有用的:预先分配该大小,如果需要,则采用摊销分配。


一条备注:R不擅长循环。对于小型循环,例如在数据框中的变量或目录中的文件上,通常没有问题。但是,如果您有一个需要使用许多循环实现的长时间计算,并且无法矢量化,则R可能不是正确的工具。在某些情况下,使用C、C++、Fortran或Java编写函数可能会有所帮助:构建插件或使用Rcpp非常容易,而性能提升也很大。


谢谢!但是那并没有完全回答我的问题。你提供的例子都是在计算之前就已经知道大小了。如果你不知道需要1e4个元素怎么办? - eduardokapp
1
设置长度比附加空向量稍微快一点。即用 k <- k + k; length(a) <- k 替换 a <- c(a, numeric(k)); k <- k+k 可以获得小幅加速。 - user2554330
@eduardokapp 在最后一个例子中,你不需要事先知道,这只是为了基准测试而需要。你可以使用while循环,跟踪实际写入的项目数量,并在最后将其截断到该大小,增加向量的机制是相同的。 - user13963867
@user2554330 我不知道这个技巧,谢谢! - user13963867

0

您可以设置向量的长度,然后对特定元素进行赋值。您的代码应该如下所示:

algo2 <- function(x) {
  if (x == 1)
    return(1)

  output <- x
  index <- 1
  while(x != 1) {
    x <- compute(x)
    index <- index + 1
    if (index > length(output))
      length(output) <- 2*length(output)
    output[index] <- x
  }

  return(output[seq_len(index)])
}

这会有所不同,虽然在您的示例中并不是很大,因为所有对compute()(和return()!)的调用都相当昂贵。如果将该计算折叠到algo中,您将看到更多的改进。您还可以将output初始化为一个长度,这个长度可能对于大多数情况来说已经足够好了,并且很少需要加倍。


没错,我的例子只是为了叙述方便而已。所以你的想法是采用传统的“需要更大的数组时将数组扩大一倍”的策略,对吧? - eduardokapp
是的,尽管我想这取决于问题,但那可能是增加的最佳方式。但主要要点是限制增加的数量,使用 length( ) <- 来设置它们,并分配到特定位置。 - user2554330

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接