2

我有一个相当大的向量（长度> 500,000）。它包含一堆NA穿插的，1并且总是保证它以 . 开头1。

我想根据对另一个向量（与 长度相同）的连续索引的比较操作，NA用v1替换其中的一些。1v2v1

是否有一种有效的方法可以在矢量化符号中执行此操作，以便在低级实现中完成循环？也许使用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
4

4 回答 4

3

可能有一个更快的解决方案，但这是我在几分钟内能想到的最好的解决方案。对于小向量，我的解决方案比 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
于 2012-07-09T14:28:34.097 回答
1

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

于 2012-07-09T12:18:56.203 回答
1

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

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]
}
于 2012-07-09T11:39:54.877 回答
1

使用编译器包可以大大加快函数 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 慢。

于 2012-07-10T09:59:49.610 回答