R中的子集xts时间序列对象答案

【问题标题】：Subset xts time-series object in RR中的子集xts时间序列对象
【发布时间】：2023-05-03 19:51:01
【问题描述】：

我有像这样的某些月份的时间序列 xts 对象

library(xts)
  seq<- seq(as.POSIXct("2015-09-01"),as.POSIXct("2015-09-04"), by = "30 mins")
  ob<- xts(data.frame(power=1:(length(seq))),seq)

现在，对应于每个观察（比如A），我想计算最后两个小时观察的平均值。因此，对应于每个观察（A）我需要计算两个小时前发生的观察的索引到A，比如说它是B。然后我可以计算A 和B 之间观察值的平均值。因此，

i=10 # dummy
ind_cur<- index(ob[i,]) # index of current observation
ind_back <- ind_cur - 3600 * 2 # index of 2 hours back observation

有了这些索引，我将 ob 子集化为

 ob['ind_cur/ind_back']

这会导致以下错误：

Error in if (length(c(year, month, day, hour, min, sec)) == 6 && c(year,  : 
  missing value where TRUE/FALSE needed
In addition: Warning messages:
1: In as_numeric(YYYY) : NAs introduced by coercion
2: In as_numeric(MM) : NAs introduced by coercion
3: In as_numeric(DD) : NAs introduced by coercion
4: In as_numeric(YYYY) : NAs introduced by coercion
5: In as_numeric(MM) : NAs introduced by coercion
6: In as_numeric(DD) : NAs introduced by coercion

谁能帮我把ob子集化！在link找到了一个相关的问题，但是还不足以解决这个问题。

更新预期输出显示为

2015-09-01 00:00:00     1   NA # as I don't have previous data
2015-09-01 00:30:00     2   NA
2015-09-01 01:00:00     3   NA
2015-09-01 01:30:00     4   NA
2015-09-01 02:00:00     5   10/4 # mean of prevous 4 observations (last two hours)
2015-09-01 02:30:00     6   14/4  
2015-09-01 03:00:00     7   18/4

【问题讨论】：

预期输出是什么？
这和移动平均线不一样吗？为此，您可以使用 TTR 包中定义的 SMA 函数。此外，在您当前的实现中，索引 'ind_cur/ind_back' 将被视为字符串文字，不会扩展到实际日期。

标签： r time-series xts

【解决方案1】：

这是一个一般难以解决的问题，因此您需要推出自己的解决方案。最简单的方法是使用 window 通过重叠 2 小时间隔进行子集化。

# initialize a result object
ob2 <- ob * NA_real_
# loop over all rows and calculate 2-hour mean
for(i in 2:nrow(ob)) {
  ix <- index(ob)[i]
  ob2[i] <- mean(window(ob, start=ix-3600*2, end=ix))
}
# set incomplete 2-hour intervals to NA
is.na(ob2) <- which(index(ob2) < start(ob2)+3600*2)

【讨论】：

【解决方案2】：

我们可以将rollapply() 包与lag() 结合使用，以将生成的滚动mean 偏移一行。

rollapply(lag(ob), 4, mean)
#                    power
#2015-09-01 00:00:00    NA
#2015-09-01 00:30:00    NA
#2015-09-01 01:00:00    NA
#2015-09-01 01:30:00    NA
#2015-09-01 02:00:00   2.5
#2015-09-01 02:30:00   3.5
#2015-09-01 03:00:00   4.5

# Or if you want it as new variable in your xts object
ob$mean <- rollapply(lag(ob),4,mean)

【讨论】：

谢谢。它解决了这个目的，但我仍然想使用xts 索引来解决这个问题，因为有时我的数据不是连续的。缺少一些读数。使用这种方法，我将无法找出计算错误。
请注意，上面的代码调用的是 rollapply.xts 而不是 rollapply.zoo。也不清楚关于连续的评论是什么意思——可能是一种误解。
根据定义，您不能将zooreg 用于不规则的时间序列。 @G.Grothendieck：我认为他们的意思是他们的实际系列不规则，所以 2 小时的间隔不均匀。
@HaroonRashid：如果缺少某些读数，您可以将对象与具有常规索引的空 xts 对象合并，以添加所需的缺失值。

【解决方案3】：

基于对“预期输出”问题的更新和 R.S. 的评论：

library(TTR)
head(SMA(ob$power, 4))  # 2 hour moving average

结果

                    SMA
2015-09-01 00:00:00  NA
2015-09-01 00:30:00  NA
2015-09-01 01:00:00  NA
2015-09-01 01:30:00 2.5
2015-09-01 02:00:00 3.5
2015-09-01 02:30:00 4.5

这假设所讨论的时间间隔为 30 分钟。

看起来更像预期输出：

lag(head(SMA(ob$power, 4),7))

                    SMA
2015-09-01 00:00:00  NA
2015-09-01 00:30:00  NA
2015-09-01 01:00:00  NA
2015-09-01 01:30:00  NA
2015-09-01 02:00:00 2.5
2015-09-01 02:30:00 3.5
2015-09-01 03:00:00 4.5

【讨论】：

【解决方案4】：

包data.table 提供滚动功能，适用于单个时间序列和多个时间序列：

head(

    as.data.table(ob)[, roll_power := frollmean(power, 4, align = 'right')]
)

# at the end of a 4 1/2 hour lag

                 index power roll_power
1: 2015-09-01 00:00:00     1         NA
2: 2015-09-01 00:30:00     2         NA
3: 2015-09-01 01:00:00     3         NA
4: 2015-09-01 01:30:00     4        2.5 # the rolling mean covers this, and preceding rows
5: 2015-09-01 02:00:00     5        3.5
6: 2015-09-01 02:30:00     6        4.5

【讨论】：