如何有条件地检查和替换 xts 对象中的数据？答案

【问题标题】：How to conditionally check and replace data in xts object?如何有条件地检查和替换 xts 对象中的数据？
【发布时间】：2020-04-16 19:33:11
【问题描述】：

这是一个可重现的数据集。问题是在一系列 NA 之间找到 1 或 2 个连续的非 NA 值并将它们分配为 NA。如果超过 2 个，则无需执行任何操作。

set.seed(55)
data <- rnorm(10)
dates <- as.POSIXct("2019-03-18 10:30:00", tz = "CET") + 0:9*60

R <- xts(x = data, order.by = dates)
colnames(R) <- "R-factor"
R[c(1, 3, 6, 10)] <- NA
R

输出：

                        R-factor
2019-03-18 10:30:00           NA
2019-03-18 10:31:00 -1.812376850
2019-03-18 10:32:00           NA
2019-03-18 10:33:00 -1.119221005
2019-03-18 10:34:00  0.001908206
2019-03-18 10:35:00           NA
2019-03-18 10:36:00 -0.505343855
2019-03-18 10:37:00 -0.099234393
2019-03-18 10:38:00  0.305353199
2019-03-18 10:39:00           NA

预期结果：

                        R-factor
2019-03-18 10:30:00           NA
2019-03-18 10:31:00           NA
2019-03-18 10:32:00           NA
2019-03-18 10:33:00           NA
2019-03-18 10:34:00           NA
2019-03-18 10:35:00           NA
2019-03-18 10:36:00 -0.505343855
2019-03-18 10:37:00 -0.099234393
2019-03-18 10:38:00  0.305353199
2019-03-18 10:39:00           NA

我编写了一个带有 for 循环的函数，该函数适用于小型数据集，但速度极慢。原始数据由 100,000+ 个数据点组成，超过 10 分钟后该函数无法执行

任何人都可以帮助我避免循环以使其更快吗？

【问题讨论】：

标签： r for-loop timestamp time-series xts

【解决方案1】：

创建一个函数Fillin，如果长度小于或等于 3，则返回 NA（如果第一个元素不是 NA，则返回 2，以便我们可以处理第一个组，即使它不以 NA 开头）否则返回其参数。使用cumsum 对运行进行分组并将Fillin 应用于每个组。

Fillin <- function(x) if (length(x) <= 3 - !is.na(x[1])) NA else x
Rc <- coredata(R)
R[] <- ave(Rc, cumsum(is.na(Rc)), FUN = Fillin)

给予：

> R
                       R-factor
2019-03-18 10:30:00          NA
2019-03-18 10:31:00          NA
2019-03-18 10:32:00          NA
2019-03-18 10:33:00          NA
2019-03-18 10:34:00          NA
2019-03-18 10:35:00          NA
2019-03-18 10:36:00 -0.50534386
2019-03-18 10:37:00 -0.09923439
2019-03-18 10:38:00  0.30535320
2019-03-18 10:39:00          NA

性能

此解决方案的运行速度与使用 rle 的解决方案大致相同。

library(microbenchmark)

microbenchmark(
  Fill = { Fillin <- function(x) if (length(x) <= 3 - !is.na(x[1])) NA else x
    Rc <- coredata(R)
    R[] <- ave(Rc, cumsum(is.na(Rc)), FUN = Fillin)
  },
  RLrep = { rleR <-  rle(c(is.na(R[,1]))) 
    is.na(R) <- with(rleR,  rep(lengths < 3 , lengths ) )
  }
)

给予：

Unit: microseconds
  expr   min    lq    mean median     uq    max neval cld
  Fill 490.9 509.5 626.550  527.7 596.45 3411.1   100   a
 RLrep 523.5 540.8 604.061  550.8 592.00 1244.4   100   a

【讨论】：

【解决方案2】：

也许可以根据Distance from the closest non NA value in a dataframe试试这个

library(tidyverse)

set.seed(55)
x <- 100000
data <- rnorm(x)
dates <- as.POSIXct("2019-03-18 10:30:00", tz = "CET") + (seq_len(x))*60
time_table1 <- tibble(time = dates,data = data)
time_table <- time_table1 %>% 
  mutate(random = rnorm(x),
         new = if_else(random > data,NA_real_,data)) %>% 
  select(-data,-random) %>% 
  rename(data= new)



lengths_na <- time_table$data %>% is.na %>% rle  %>% pluck('lengths')

the_operation <- . %>% 
  mutate(lengths_na =lengths_na %>% seq_along %>% rep(lengths_na)) %>% 
  group_by(lengths_na) %>%
  add_tally() %>%
  ungroup() %>% 
  mutate(replace_sequence = if_else(condition = n < 3,true = NA_real_,false = data))

microbenchmark::microbenchmark(time_table %>% the_operation)

效果还不错

Unit: milliseconds
                         expr      min       lq     mean  median       uq      max neval
 time_table %>% the_operation 141.9009 176.2988 203.3744 190.183 214.1691 412.3161   100

也许这更容易阅读

library(tidyverse)

set.seed(55)

# Create the data

x <- 100
data <- rnorm(x)
dates <- as.POSIXct("2019-03-18 10:30:00", tz = "CET") + (seq_len(x))*60
time_table1 <- tibble(time = dates,data = data)

# Fake some na's
time_table <- time_table1 %>% 
  mutate(random = rnorm(x),
         new = if_else(random > data,NA_real_,data)) %>%
  select(-data,-random) %>% 
  rename(data= new)


# The rle function counts the occurrences of the same value in a vector,
# We create a T/F vector using is.na function
# meaning that we can count the lenght of sequences with or without na's
lengths_na <- time_table$data %>% is.na %>% rle  %>% pluck('lengths')

# This operation here can be done outside of the df
new_col <- lengths_na %>%
  seq_along %>% # Counts to the size of this vector
  rep(lengths_na) # Reps the lengths of the sequences populating the vector

result <- time_table %>%
  mutate(new_col =new_col) %>% 
  group_by(new_col) %>% # Operates the logic on this group look into the tidyverse
  add_tally() %>% # Counts how many instance there are on each group 
  ungroup() %>% # Not actually needed but good manners
  mutate(replace_sequence = if_else(condition = n < 3,true = NA_real_,false = data))

【讨论】：

@Bruno，它几乎可以按要求工作。几乎是因为，它还取代了一系列非 NA 之间的单独或一对 NA（这对我来说可能不是问题）。作为一个初学者，这对我来说理解起来有点复杂，但它实际上在原始大型数据集上运行得非常快。所以谢谢你
嗨@user11841472，您可以使 if_else 语句更复杂以忽略 na，但我认为这不会改变最终结果，如果您在 rstudio 上，我也会处理一些 cmets滥用 f1 键阅读您不知道的功能，也许更改您的 SO 名称？ user11841472 感觉没有人情味
不断收到"Error in rename(., data = new) : unused argument (data = new)"
我可以在我现有的工作空间中复制它，但重新启动使其成功。但是它不适用于 xts 对象。
我使用 fortify.zoo() 将 xts 中的原始数据转换为 df 并且它工作了 @42

【解决方案3】：

我想，还有更优雅的解决方案，但这会将时间缩短一半

    R_df=as.data.frame(R)

    R_df$shift_1=c(R_df$`R-factor`[-1],NA) #shift value one up
    R_df$shift_2=c(NA,R_df$`R-factor`[-nrow(R_df)]) #shift value one down

# create new filtered variable
    R_df$`R-factor_new`=ifelse(is.na(R_df$`R-factor`),NA,
                               ifelse((!is.na(R_df$shift_1))|(!is.na(R_df$shift_2)),
                                      R_df$`R-factor`,NA)

>                 test replications elapsed relative user.self sys.self user.child sys.child
>     2 ifelseapproach         1000    0.83    1.000      0.65     0.19         NA        NA
>     1       original         1000    1.81    2.181      1.76     0.01         NA        NA

【讨论】：

您好，此解决方案仅替换一个非 NA，但遗憾的是缺少替换 2 个连续的非 NA 值。
好吧，你在编辑之前在原始问题中发布的函数removeErrorFun <- function(temp){ for (i in 1:length(temp)-1){ if (is.na(temp[i-1]$R-factor) && is.na(temp[i+1]$R-factor)){ temp[i]$R-factor` = NA } } return(temp) }` left中的两个连续值。所以我似乎误读了您的原始帖子。您可以轻松扩展我的方法，但现在您还有其他选择 ;-)
你说得对，把这个功能放在首位是我的错误。我一意识到就删除了这个功能。带来不便敬请谅解。 @TobiO

【解决方案4】：

这可能比提供的大多数其他解决方案更快。 rep 函数本质上是 rle 函数的逆函数。它采用两个向量参数，并将第一个的值的计数扩展为第二个的长度，这允许基于运行长度进行测试，然后替换为is.na <-。实际上有两个不同的函数： rle(x) 它返回一个长度为 (x) 的逻辑向量，然后是 is.na(x)<- 它根据向量右侧的逻辑值将 NA 分配给 x 中的项目那个函数。：

rleR <- rle(c(is.na(R[,1]))) #get the position and lengths of nonNA's and NA's
is.na(R) <- with(rleR,  rep(lengths < 3 , lengths ) ) #set NAs
#--------------
> R
                       R-factor
2019-03-18 10:30:00          NA
2019-03-18 10:31:00          NA
2019-03-18 10:32:00          NA
2019-03-18 10:33:00          NA
2019-03-18 10:34:00          NA
2019-03-18 10:35:00          NA
2019-03-18 10:36:00 -0.50534386
2019-03-18 10:37:00 -0.09923439
2019-03-18 10:38:00  0.30535320
2019-03-18 10:39:00          NA
Warning message:
timezone of object (CET) is different than current timezone (). 


microbenchmark(
 Fill = {Fillin <- function(x) if (length(x) <= 3 - !is.na(x[1])) NA else x
ave(R, cumsum(is.na(R)), FUN = Fillin)}, 
 RLrep = {rleR <-  rle(c(is.na(R[,1])))
     is.na(R) <- with(rleR,  rep(lengths < 3 , lengths ) )})
#----------------------
Unit: microseconds
  expr      min        lq      mean    median        uq      max neval cld
  Fill 1668.788 1784.6275 1942.5261 1844.5825 2005.0960 4911.762   100   b
 RLrep  102.174  113.9565  144.3477  131.4735  156.6715  368.665   100  a

【讨论】：

我稍微修改了我的代码。请参阅我的帖子末尾的基准。现在它们的运行速度大致相同，其中一个具有更快的平均时间，另一个具有更快的平均时间。
嗯，我认为“我的”似乎更加矢量化，但我很确定我多年前在 Rhelp 上从你那里学到了这个策略，所以我不应该太自大。
@42- 这段代码工作正常，但我想我是完全理解你的代码的新手。只是为了了解你做了什么，你能详细说明一下这条线吗？ is.na(R)
你需要意识到is.na<-实际上是一个选择性赋值函数，而不是is.na和<-的组合。请参阅上面的编辑。