【问题标题】:R fill vector efficientlyR 有效地填充向量
【发布时间】:2023-03-03 18:06:01
【问题描述】:

我有一个相当大的向量(长度>500,000)。它包含一堆NA 穿插1 并且始终保证它以1 开头。

我想将v1 中的一些NA 替换为1,基于对另一个向量v2(与v1 长度相同)的连续索引的比较操作。

有没有一种有效的方法以矢量化表示法执行此操作,以便在低级实现中完成循环?也许使用ifelse

下面的可重现示例:

v1<-c(1,NA,NA,NA,1,NA,NA,NA,NA,NA,1,NA,NA,1,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,1)
v2<-c(10,10,10,9,10,9,9,9,9,9,10,10,10,11,8,12,12,12,12,12,12,12,12,12,12,13)
# goal is to fill through v1 in such a way that whenever 
# v1[i] == NA and v1[i-1] == 1 and v2[i] == v2[i-1], then v1[i] == 1
MM<-data.frame(v1,v2)
for (i in 2:length(v1)){ 
    # conditions: v1[i-1] == 1; v1[i]==NA; v2[i]==v2[i-1]
    if (!is.na(v1[i-1]) && is.na(v1[i]) && v2[i]==v2[i-1]){
        v1[i]<-1
    }
}
MM$v1_altered<-v1
MM

【问题讨论】:

  • 你能提供一个 v2 的例子吗? IE。 reproducible example...
  • @JoshuaUlrich 我编辑了我的原始帖子,添加了一个可重复的示例。复制粘贴应该可以,谢谢
  • 您的可重现示例与您在 v1v2 上运行的初始、不可重现的示例不同。哪个包含你想要的输出?
  • @JoshuaUlrich 可重现的例子是正在考虑的问题,很抱歉造成混淆——最初的不可重现是问题的本质,但不是精确的规范
  • R populating a vector的可能重复

标签: r loops


【解决方案1】:

可能有一个更快的解决方案,但这是我在几分钟内能想到的最佳解决方案。对于小向量,我的解决方案比 OP 慢,但对于更大的向量,我的解决方案却越来越快。

library(zoo)  # for na.locf
library(rbenchmark)

v1<-c(1,NA,NA,NA,1,NA,NA,NA,NA,NA,1,NA,NA,1,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,1)
v2<-c(10,10,10,9,10,9,9,9,9,9,10,10,10,11,8,12,12,12,12,12,12,12,12,12,12,13)
V1 <- rep(v1, each=20000)  # 520,000 observations
V2 <- rep(v2, each=20000)  # 520,000 observations

fun1 <- function(v1,v2) {
  for (i in 2:length(v1)){ 
    if (!is.na(v1[i-1]) && is.na(v1[i]) && v2[i]==v2[i-1]){
      v1[i]<-1
    }
  }
  v1
}
fun2 <- function(v1,v2) {
  # create groups in which we need to assess missing values
  d <- cumsum(as.logical(c(0,diff(v2))))
  # for each group, carry the first obs forward
  ave(v1, d, FUN=function(x) na.locf(x, na.rm=FALSE))
}
all.equal(fun1(V1,V2), fun2(V1,V2))
# [1] TRUE
benchmark(fun1(V1,V2), fun2(V1,V2))
#           test replications elapsed relative user.self sys.self
# 1 fun1(V1, V2)          100  194.29 6.113593    192.72     0.17
# 2 fun2(V1, V2)          100   31.78 1.000000     30.74     0.95

【讨论】:

    【解决方案2】:

    矢量化解决方案如下所示:

    v1[-1] <- ifelse(diff(v2), 0, v1[-length(v1)])
    

    但上述方法不起作用,我认为您无法避免显式循环,因为如果我理解正确,您想传播新值。那么,怎么样:

    cmp <- diff(v2)
    for (i in 2:length(v1)){
        v1[i] <- if(cmp[i-1]) 0 else v1[i-1]
    }
    

    【讨论】:

      【解决方案3】:

      它可能不会更快,但v1[i] &lt;- v1[i-1] * (cmp[i-1] == 0) 避免了所有显式的“if”调用。我现在无法对其进行测试,但您可以尝试@James 解决方案与循环遍历此表单,例如长度为 1e4 的向量,看看哪个执行得更快。

      【讨论】:

        【解决方案4】:

        使用编译器包可以大大加快函数 fun1 的速度。 使用 Joshua 提供的代码并使用编译器包对其进行扩展:

        library(zoo)  # for na.locf
        library(rbenchmark)
        library(compiler)
        
        v1 <- c(1,NA,NA,NA,1,NA,NA,NA,NA,NA,1,NA,NA,1,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,1)
        v2 <- c(10,10,10,9,10,9,9,9,9,9,10,10,10,11,8,12,12,12,12,12,12,12,12,12,12,13)
        
        fun1 <- function(v1,v2) {
            for (i in 2:length(v1)){
                if (!is.na(v1[i-1]) && is.na(v1[i]) && v2[i]==v2[i-1]){
                    v1[i]<-1
                }
            }
            v1
        }
        
        fun2 <- function(v1,v2) {
            # create groups in which we need to assess missing values
            d <- cumsum(as.logical(c(0,diff(v2))))
            # for each group, carry the first obs forward
            ave(v1, d, FUN=function(x) na.locf(x, na.rm=FALSE))
        }
        
        fun3 <- cmpfun(fun1)
        
        fun1(v1,v2)
        fun2(v1,v2)
        all.equal(fun1(v1,v2), fun2(v1,v2))
        all.equal(fun1(v1,v2), fun3(v1,v2))
        
        Nrep <- 1000
        
        V1 <- rep(v1, each=Nrep)
        V2 <- rep(v2, each=Nrep)
        all.equal(fun1(V1,V2), fun2(V1,V2))
        all.equal(fun1(V1,V2), fun3(V1,V2))
        
        benchmark(fun1(V1,V2), fun2(V1,V2), fun3(V1,V2))
        

        我们得到以下结果

        benchmark(fun1(V1,V2), fun2(V1,V2), fun3(V1,V2))
                  test replications elapsed relative user.self sys.self user.child
        1 fun1(V1, V2)          100  12.252 5.706567    12.190    0.045          0
        2 fun2(V1, V2)          100   2.147 1.000000     2.133    0.013          0
        3 fun3(V1, V2)          100   3.702 1.724266     3.644    0.023          0
        

        所以编译后的 fun1 比原来的 fun1 快很多,但仍然比 fun2 慢。

        【讨论】:

        • fun3fun2 的编译版本,而不是 fun1。在我的示例中,fun1 的编译版本仍然比 fun2 慢约 2 倍。
        • @Joshua。愚蠢的错误。 fun1 的编译版本确实几乎是 fun2 的两倍,比原来的 fun1 快三倍多。
        • 我已经更正了答案,现在它显示了预期的内容。
        猜你喜欢
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 2018-11-28
        相关资源
        最近更新 更多