随窗口变化的移动平均线答案

【问题标题】：Moving average with changing window随窗口变化的移动平均线
【发布时间】：2018-08-09 05:48:42
【问题描述】：

我想计算 R 中一个变量的移动平均值，窗口大小会发生变化。更具体地说：移动平均值应在三年内计算，但数据（时间序列）的可用频率更高，并且每个三年窗口的窗口大小可能会有所不同。

假设以下数据集：

library(data.table)
set.seed(1)   # reproduceable data
dataset <- data.table(ID=c(rep("A",2208),rep("B",2208)),
x = c(rnorm(2208*2)), time=c(seq(as.Date("1988/03/15"),
as.Date("2000/04/16"), "day"),seq(as.Date("1988/03/15"),
as.Date("2000/04/16"), "day")))

应为两个 ID（人 A 和 B）计算变量 x 的三年移动平均值。这可以用zoo 和datatable 来完成吗？但是任何解决方案都可以。

请注意，我知道如何使用固定的窗口大小来执行此操作，这里的问题是窗口大小不同。

【问题讨论】：

您所说的每个三年窗口的窗口大小可能会有所不同是什么意思？
我认为 OP 意味着窗口将跨越三年，但这些年内的观察次数是可变的。
如果数据集大小不是问题，则创建一个 data.table，从 min date 开始，到 max date 并加入您的数据集，然后您可以使用固定的滑动窗口大小。
@Lyngbakr 这就是我的意思，谢谢。不幸的是，数据集大小是一个问题，但即使不是，窗口大小不规则的问题也会存在，因为在给定的一年中不一定有与另一年相同的天数（即在闰年总是有与正常年份不同）？
我不知道它是否有效，但是 dataset %>% group_by(ID) %>% mutate(avg=rollapplyr(x,mean,width=3*365,fill=NA)) 怎么样，因为您已经有没有缺失值的每日测量值。

标签： r

【解决方案1】：

如果我理解正确的话，OP 想要跨越 3 年。由于这可能包括闰年，因此窗口大小可以是 1095 天或 1096 天。

这可以通过在非等值连接中聚合以及lubridate的回滚日期算法来解决。

library(data.table)
library(lubridate)
# create 3 years windows for each ID for later non-equi join
win <- dataset[, CJ(ID = ID, start = time, unique = TRUE)][
  # make sure to pick
  , end := start %m+% years(3) - days(1)][
    # remove windows which end out of date range
    end <= max(start)]
win

      ID      start        end
   1:  A 1988-03-15 1991-03-14
   2:  A 1988-03-16 1991-03-15
   3:  A 1988-03-17 1991-03-16
   4:  A 1988-03-18 1991-03-17
   5:  A 1988-03-19 1991-03-18
  ---                         
6638:  B 1997-04-13 2000-04-12
6639:  B 1997-04-14 2000-04-13
6640:  B 1997-04-15 2000-04-14
6641:  B 1997-04-16 2000-04-15
6642:  B 1997-04-17 2000-04-16

# check window lengths
win[, .N, by = .(days = end - start + 1L)]

        days    N
1: 1095 days 2166
2: 1096 days 4476

# see what happens in leap years
win[leap_year(start) & month(start) == 2 & day(start) %in% 28:29, 
  .(start, end, days = end - start + 1L)]

        start        end      days
1: 1992-02-28 1995-02-27 1096 days
2: 1992-02-29 1995-02-27 1095 days
3: 1996-02-28 1999-02-27 1096 days
4: 1996-02-29 1999-02-27 1095 days
5: 1992-02-28 1995-02-27 1096 days
6: 1992-02-29 1995-02-27 1095 days
7: 1996-02-28 1999-02-27 1096 days
8: 1996-02-29 1999-02-27 1095 days

win[leap_year(end) & month(end) == 2 & day(end) %in% 28:29,
    .(start, end, days = end - start + 1L)]

        start        end      days
1: 1989-03-01 1992-02-29 1096 days
2: 1993-03-01 1996-02-29 1096 days
3: 1997-03-01 2000-02-29 1096 days
4: 1989-03-01 1992-02-29 1096 days
5: 1993-03-01 1996-02-29 1096 days
6: 1997-03-01 2000-02-29 1096 days

# aggregate in a non-equi-join
dataset[win, on = .(ID, time >= start, time <= end), by = .EACHI, .(avg = mean(x))]

      ID       time       time         avg
   1:  A 1988-03-15 1991-03-14 -0.01184078
   2:  A 1988-03-16 1991-03-15 -0.01317813
   3:  A 1988-03-17 1991-03-16 -0.01179571
   4:  A 1988-03-18 1991-03-17 -0.01006100
   5:  A 1988-03-19 1991-03-18 -0.01221798
  ---                                     
6638:  B 1997-04-13 2000-04-12 -0.03412214
6639:  B 1997-04-14 2000-04-13 -0.03604176
6640:  B 1997-04-15 2000-04-14 -0.03556291
6641:  B 1997-04-16 2000-04-15 -0.03392185
6642:  B 1997-04-17 2000-04-16 -0.03393674

【讨论】：

感谢这对我来说很完美。只有一个问题：你能不能用 datatable 中的常用 merge 语法来表述这个？
不，像 non-equi join 和 grouping by each i 这样的好东西只能在 [.data.table 语法中使用。 data.table 的 merge() 或多或少只是基本 R 的 merge() 的更快实现。 a section in the FAQ 解释了这些差异。
我明白了，谢谢。我问是因为在我看来merge 的语法更干净。很难记住 X[Y] 和 Y[X] imho 之间的区别。但是解决方案很棒，所以我想我需要习惯这一点。
是的，merge() 更简单但功能更弱。 X[Y] 的行为类似于右连接，即获取 Y 的所有行并尝试在 X 中查找匹配项。
有包fuzzyjoin 也可以工作，并且语法接近merge（甚至更接近dplyr 连接），但我相信它会更慢。

【解决方案2】：

正如@A.Suliman 在 cmets 中提到的那样，您的示例数据具有固定的窗口宽度，但我们假设您的真实数据没有您的文本所说的那样。

rollapply 中的参数width 不必是常数，所以我们可以先计算宽度，确保窗口左对齐，然后运行rollaply。

library(zoo)
library(tidyverse)
dataset %>% 
  arrange(ID,time) %>%
  group_by(ID) %>%
  mutate(avg = rollapply(x, FUN = mean, align = "left",
                         width = map_dbl(time, ~which.max(time[time < .x + 3*365.25])) - row_number()+1))
#                          
# # A tibble: 8,832 x 4
# # Groups:   ID [2]
#         ID       x       time    avg
#     <fctr>   <dbl>     <date>  <dbl>
#   1      A -0.0258 1988-03-15 0.0109
#   2      A -0.0258 1988-03-15 0.0109
#   3      A -0.1562 1988-03-16 0.0115
#   4      A -0.1562 1988-03-16 0.0115
#   5      A  0.8193 1988-03-17 0.0115
#   6      A  0.8193 1988-03-17 0.0112
#   7      A -1.1136 1988-03-18 0.0102
#   8      A -1.1136 1988-03-18 0.0107
#   9      A -0.9336 1988-03-19 0.0105
#  10      A -0.9336 1988-03-19 0.0109

【讨论】：

我很着急，所以没有费心将语法翻译成 data.table 但我认为这应该是微不足道的。解决方案未深入检查，谨慎对待！
谢谢，很有趣。 map_dbl 似乎是一个 purrr 函数，所以我想我需要那个包。你介意解释一下~which.max(time[time < .x + 3*365.25])) - row_number()+1) 背后的想法吗？
您可以使用sapply 代替map_dbl，使用标准函数表示法。您复制的代码背后的想法是找到给定日期 + 3 年之前可用的最后日期的索引，从中减去给定日期的索引 (+1) 以获得时间窗口，并在 width 中使用此窗口论据
乍一看，我的解决方案似乎和@Uwe 的原理几乎相同。他的更干净，而且所有 data.table 都是如此，所以你最好专注于它。