根据索引有条件地替换向量的元素答案

【问题标题】：Conditionally replace elements of a vector based on an index根据索引有条件地替换向量的元素
【发布时间】：2017-12-08 01:37:44
【问题描述】：

最好用一个例子来解释。

我有一个来自data.frame 的向量或列，名为vec：

vec <- c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA)

我想要一个矢量化过程（不是for 循环）在观察到1 时更改三个尾随NA。

结束向量是：

c(NA, NA, 1, 1, 1, 1, NA, 1, 1, 1, 1, NA, NA, NA)

如果我们有：

vec <- c(NA, NA, 1, NA, 1, NA, NA, 1, NA, NA, NA, NA, NA, NA)

结束向量如下所示：

c(NA, NA, 1, 1, 1, 1, 1, 1, 1, 1, 1, NA, NA, NA)

一个写得很糟糕的解决方案是：

vec2 <- vec
for(i in index(v)){
  if(!is.na(v[i])) vec2[i] <- 1
  if(i>3){
    if(!is.na(vec[i-1])) vec2[i] <- 1
    if(!is.na(vec[i-2])) vec2[i] <- 1
    if(!is.na(vec[i-3])) vec2[i] <- 1
  }
  if(i==3){
    if(!is.na(vec[i-1])) vec2[i] <- 1
    if(!is.na(vec[i-2])) vec2[i] <- 1
  }
  if(i==2){
    if(!is.na(vec[i-1])) vec2[i] <- 1
  }
}

【问题讨论】：

Welp，R 中没有矢量化的seq。我能想到的最接近的方法是使用Vectorize，如seq2 <- Vectorize(seq.default, vectorize.args = c("from", "to")) ; indx <- which(!is.na(vec)) ; vec[seq2(indx, indx + 3)] <- 1。取自here。但话又说回来，Vectorize 只是mapply。所以你可以用 Rcpp 写一个 vecotrized seq 可能
Fill NA in a time series only to a limited number的可能重复
@raymkchow 谢谢没有找到那个帖子。让每个人都想办法解决，我除了跑得更快的那个。
@dimitris_ps 如果接受的答案是您正在寻找的（不是矢量化解决方案），那么这与我在原始评论中发布的完全一样，这只是一个骗局。跨度>

标签： r

【解决方案1】：

另一种选择：

`[<-`(vec,c(outer(which(vec==1),1:3,"+")),1)
# [1] NA NA  1  1  1  1 NA  1  1  1  1 NA NA NA

虽然上面的例子适用于上面的例子，但如果在最后一个位置找到 1，它会拉伸 vec 的长度。最好做一个简单的检查并包装成一个函数：

threeNAs<-function(vec) {
    ind<-c(outer(which(vec==1),1:3,"+"))
    ind<-ind[ind<=length(vec)]
    `[<-`(vec,ind,1)
}

【讨论】：

一个较小的替代方案如下：spots <- unique(c(outer(which(!is.na(vec)), 1:3, "+"))); vec[spots[spots <= length(vec)]] <-1。这应该考虑重复索引以及保持原始长度。

【解决方案2】：

另一种快速解决方案：

vec[rep(which(vec == 1), each = 3) + c(1:3)] <- 1

给出：

> vec
 [1] NA NA  1  1  1  1 NA  1  1  1  1 NA NA NA

基准测试仅在对较大的数据集进行时才真正有用。具有 10k 大向量的基准测试和几个已发布的解决方案：

library(microbenchmark)

microbenchmark(ans.jaap = {vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4); 
                           vec[rep(which(vec == 1), each = 3) + c(1:3)] <- 1},
               ans.989 = {vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
                          r <- which(vec==1);
                          vec[c(mapply(seq, r, r+3))] <- 1},
               ans.sotos = {vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
                            vec[unique(as.vector(t(sapply(which(vec == 1), function(i) seq(i+1, length.out = 3)))))] <- 1},
               ans.gregor = {vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
                             vec[is.na(vec)] <- 0;
                             n <- length(vec);
                             vec <- vec + c(0, vec[1:(n-1)]) + c(0, 0, vec[1:(n-2)]) + c(0, 0, 0, vec[1:(n-3)]);
                             vec[vec == 0] <- NA},
               ans.moody = {vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
                            output <- sapply(1:length(vec),function(i){any(!is.na(vec[max(0,i-3):i]))});
                            output[output] <- 1;
                            output[output==0] <- NA},
               ans.nicola = {vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
                             `[<-`(vec,c(outer(which(vec==1),1:3,"+")),1)})

它给出了以下基准：

Unit: microseconds
       expr        min         lq       mean     median         uq        max neval   cld
   ans.jaap   1778.905   1937.414   3064.686   2100.595   2257.695  86233.593   100 a    
    ans.989  87688.166  89638.133  96992.231  90986.269  93326.393 182431.366   100   c  
  ans.sotos 125344.157 127968.113 132386.664 130117.438 132951.380 214460.174   100    d 
 ans.gregor   4036.642   5824.474  10861.373   6533.791   7654.587  87806.955   100  b   
  ans.moody 173146.810 178369.220 183698.670 180318.799 184000.062 264892.878   100     e
 ans.nicola    966.927   1390.486   1723.395   1604.037   1904.695   3310.203   100 a

【讨论】：

感谢您将我纳入基准测试并取消删除（不确定您是您还是@DavidArenburg：如果是的话，谢谢大卫）我的回答。我删除了它，因为which() 部分确实很常见。接受的答案的其余部分实施很差，我发表了一条评论，希望进行编辑，这会使我的答案毫无用处。

【解决方案3】：

如果不是用 C 语言编写的循环，那么真正的“向量化”是什么？

这是一个基准测试良好的 C++ 循环。

vec <- c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA)

library(Rcpp)

cppFunction('NumericVector fixVec(NumericVector myVec){

    int n = myVec.size();
    int foundCount = 0;

    for(int i = 0; i < n; i++){
      if(myVec[i] == 1) foundCount = 1; 

      if(ISNA(myVec[i])){
        if(foundCount >= 1 & foundCount <= 3){
          myVec[i] = 1;
          foundCount++;
        }
      }
    }
    return myVec;
    }')

fixVec(vec)
# [1] NA NA  1  1  1  1 NA  1  1  1  1 NA NA NA

基准测试

library(microbenchmark)

microbenchmark(
      ans.jaap = {
        vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4); 
      vec[rep(which(vec == 1), each = 4) + c(0:3)] <- 1
},

    ans.nicola = {
        vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
      `[<-`(vec,c(outer(which(vec==1),0:3,"+")),1)
        },

    ans.symbolix = {
        vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4);
      vec <- fixVec(vec)
        }
)

# Unit: microseconds
# expr              min       lq      mean   median        uq       max neval
# ans.jaap     2017.789 2264.318 2905.2437 2579.315 3588.4850  4667.249   100
# ans.nicola   1242.002 1626.704 3839.4768 2095.311 3066.4795 81299.962   100
# ans.symbolix  504.577  533.426  838.5661  718.275  966.9245  2354.373   100


vec <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4)
vec <- fixVec(vec)

vec2 <- rep(c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA),1e4)
vec2[rep(which(vec2 == 1), each = 4) + c(0:3)] <- 1

identical(vec, vec2)
# [1] TRUE

【讨论】：

编写循环的语言不是那关键，主要问题是在每次迭代中评估一个R函数，而不是做一切都在 C 中。

【解决方案4】：

以下代码可以满足您的要求。它涉及“移位”向量，然后添加移位的版本

vec[is.na(vec)] <- 0                                 
n <- length(vec)                                     
vec <- vec + c(0, vec[1:(n-1)]) + c(0, 0, vec[1:(n-2)]) + c(0, 0, 0, vec[1:(n-3)])  
vec[vec == 0] <- NA                                    
vec[vec != 0] <- 1                                     

# vec                    |   0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0 ,0, 0
# c(0, vec[1:(n-1)])     | + 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0 ,0, 0
# c(0, 0, vec[1:(n-2)])  | + 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0 ,0
# c(0,0,0,vec[1:(n-3)])  | + 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0 
#                        |-------------------------------------------
#                        |   0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 0

【讨论】：

【解决方案5】：

非矢量化解决方案，但仍然是使用基 R 的另一种选择，

vec[unique(as.vector(t(sapply(which(vec == 1), function(i) seq(i+1, length.out = 3)))))] <- 1
vec
#[1] NA NA  1  1  1  1 NA  1  1  1  1 NA NA NA

vec1[unique(as.vector(t(sapply(which(vec1 == 1), function(i) seq(i+1, length.out = 3)))))] <- 1
vec1
#[1] NA NA  1  1  1  1  1  1  1  1  1 NA NA NA

【讨论】：

【解决方案6】：

这个怎么样：

r <- which(vec==1)
vec[c(mapply(seq, r, r+3))] <- 1

例子：

vec <- c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, NA)
#[1] NA NA  1  1  1  1 NA  1  1  1  1 NA NA NA

vec <- c(NA, NA, 1, NA, 1, NA, NA, 1, NA, NA, NA, NA, NA, NA)
#[1] NA NA  1  1  1  1  1  1  1  1  1 NA NA NA

【讨论】：

这个代码将改变向量的长度，以防最后有一个1。例如vec <- c(NA, NA, 1, NA, NA, NA, NA, 1, NA, NA, NA, NA, NA, 1)
我看到你躲在标志后面，所以又来了：outer(r, 0:3, '+') 与 mapply(seq, r, r+3) 完全不同，因为 A：outer 是完全矢量化的，而 mapply 只是一个循环。 B：运行+ 比seq 便宜得多（在计算上），而且总体上是一个完全不同的概念。 @nicolas 的想法更好，满足了 vectoization 要求，甚至与您的解决方案不相似。另外，请查看基准和 ~X60 时间速度差异。
不要在回答中包含对投票行为的投诉。
@Jaap 非常正确，尽管我理解他的沮丧。它不是严格意义上的矢量化，但它是一个有效的答案，不值得那么多反对票（准确地说是 8 个）。我不明白为什么有人会投票删除它。这只是刻薄，没有任何帮助。
@JorisMeys 如果问题要求“矢量化”（我们可能会阅读快速或高效）解决方案并且有人用非矢量化解决方案回答，那是一个好（或使用 SO术语“有用”）答案？如果您认为这是一个有用的答案，请投票；那些不认为它有用的人投反对票。我看不出“有效”与此有什么关系？

【解决方案7】：

使用sapply、any 和is.na：

output <- sapply(1:length(vec),function(i){any(!is.na(vec[max(0,i-3):i]))})
output[output] <- 1
output[output==0] <- NA

【讨论】：