面板数据的滚动平均值（有一些细节）答案

【问题标题】：Rolling average for panel data (with a few details)面板数据的滚动平均值（有一些细节）
【发布时间】：2019-08-26 05:41:54
【问题描述】：

我想出了一些代码来计算面板数据的滚动平均值（数据中的一行包含一天中一个主题的值）。由于我有一些更具体的要求，代码变得相当复杂。在我看来，对于一个不太罕见的应用程序来说太复杂了。

这是我需要的：

滚动平均值（(a) 前 3 天不包括“当前”日的值的平均值，(b) 仅在有此窗口中至少有 2 个非缺失值）
尊重面板结构

不会太复杂吧？

对于 1. 我决定使用 rollapplyr() 和 mean( , na.rm = T) 来排除当天 (a) 我决定使用自制的滞后函数和 (b) if 语句。对于 2。我将所有内容都包裹在 tapply()（与 unlist()）中，以尊重面板结构。

代码示例如下：

library(zoo)

# example data (with missings)
set.seed(1)
df = data.frame(subject = rep(c("a", "b"), each = 10), day = rep(1:10, 2), value = rnorm(20))
df$value[15:17] = NA

# lag function (sensitive to "single day" subjects)
lag <- function(x, l = 1) { 
  if (length(x) > 1) (c(rep(NA, l), x[1:(length(x)-l)])) else (NA) 
} 

# calculate rolling mean
df$roll_mean3 = unlist(tapply(df$value, df$subject, 
                              FUN = function(x) lag(rollapplyr(x, width = 3, fill = NA, partial = T,
                                                               FUN = function(x) ifelse(sum(!is.na(x)) > 1, mean(x, na.rm = T), NA)))))
df

正如我所说，对于我认为并不遥远的情况，这种解决方案似乎过于复杂。

您对如何以更简单（不易出错）的方式执行此操作有什么建议吗？我是否错过了一些可以更轻松地处理面板数据的基本功能？

为了说明，我的代码的输出是：

   subject day      value   roll_mean3
1        a   1 -0.6264538           NA
2        a   2  0.1836433           NA
3        a   3 -0.8356286 -0.221405243
4        a   4  1.5952808 -0.426146366
5        a   5  0.3295078  0.314431838
6        a   6 -0.8204684  0.363053321
7        a   7  0.4874291  0.368106730
8        a   8  0.7383247 -0.001177187
9        a   9  0.5757814  0.135095124
10       a  10 -0.3053884  0.600511703
11       b   1  1.5117812           NA
12       b   2  0.3898432           NA
13       b   3 -0.6212406  0.950812202
14       b   4 -2.2146999  0.426794608
15       b   5         NA -0.815365744
16       b   6         NA -1.417970234
17       b   7         NA           NA
18       b   8  0.9438362           NA
19       b   9  0.8212212           NA
20       b  10  0.5939013  0.882528703

【问题讨论】：

您能否修改您的帖子以包含可重现的样本数据（使用set.seed 表示固定的随机种子）和相应的预期输出？这听起来像 rollapply 和 dplyr::group_by 可能会完成这项工作，但您的预期输出将有助于更好地理解您的问题。此外，dplyr 已经有一个 lag 函数（并且data.table 有一个等效的 shift 函数）。
我在问题中添加了种子语句和预期输出。
嗯。我不确定我是否理解您的预期输出。您想计算 (a) 前 3 天（不包括“当前”天）的值的 “平均值”。主题"a" 的第一个非NA 值不会从第4 行开始吗？为什么你在第 3 行有一个值？
因为一旦滚动窗口中有 2 个非缺失值，就会计算平均值。
但根据您的描述，您想计算 前 3 天（不包括当天）的 value 的平均值。这与计算平均值 “只要滚动窗口中有 2 个非缺失值” 有什么关系？似乎是两个不同的目标。

标签： r lag moving-average panel-data tapply

【解决方案1】：

使用ave 对每个主题分别运行rollapply。然后在使用rollapply 时请注意，width 可以是一个包含偏移向量（或多个向量）的列表，因此list(-seq(3)) 表示前三个元素。有关参数的更多信息，请参阅?rollapply。

Mean <- function(x) if (sum(!is.na(x)) >= 2) mean(x, na.rm = TRUE) else NA
roll <- function(x)  rollapply(x, list(-seq(3)), Mean, fill = NA, partial = TRUE)
transform(df, roll = ave(value, subject, FUN = roll))

【讨论】：

这很好。我喜欢你不必使用滞后功能。我不知道 rollapply 本身就有这样的功能（对我来说，它应该有它很有意义）。
我也喜欢你不使用 dplyr（或 data.table）。不过，这可能只是个人口味……

【解决方案2】：

关于我上面的评论，我不完全确定您期望的输出应该是什么，但也许以下是一个很好的起点：

df %>%
    group_by(subject) %>%
    mutate(roll_mean3 = rollapplyr(
        lag(value),
        width = 3,
        fill = NA,
        FUN = function(x) ifelse(sum(!is.na(x)) > 1, mean(x, na.rm = T), NA)))
## A tibble: 20 x 4
## Groups:   subject [2]
#   subject   day   value roll_mean3
#   <fct>   <int>   <dbl>      <dbl>
# 1 a           1  -0.626   NA
# 2 a           2   0.184   NA
# 3 a           3  -0.836   -0.221
# 4 a           4   1.60    -0.426
# 5 a           5   0.330    0.314
# 6 a           6  -0.820    0.363
# 7 a           7   0.487    0.368
# 8 a           8   0.738   -0.00118
# 9 a           9   0.576    0.135
#10 a          10  -0.305    0.601
#11 b           1   1.51    NA
#12 b           2   0.390   NA
#13 b           3  -0.621    0.951
#14 b           4  -2.21     0.427
#15 b           5  NA       -0.815
#16 b           6  NA       -1.42
#17 b           7  NA       NA
#18 b           8   0.944   NA
#19 b           9   0.821   NA
#20 b          10   0.594    0.883

或使用data.table

custom_mean <- function(x) ifelse(sum(!is.na(x)) > 1, mean(x, na.rm = T), NA)
setDT(df)[, roll_mean3 := rollapplyr(shift(value), width = 3, fill = NA, FUN = custom_mean), by = subject]
df
#   subject day      value   roll_mean3
#1:       a   1 -0.6264538           NA
#2:       a   2  0.1836433           NA
#3:       a   3 -0.8356286 -0.221405243
#4:       a   4  1.5952808 -0.426146366
#5:       a   5  0.3295078  0.314431838
#6:       a   6 -0.8204684  0.363053321
#7:       a   7  0.4874291  0.368106730
#8:       a   8  0.7383247 -0.001177187
#9:       a   9  0.5757814  0.135095124
#10:       a  10 -0.3053884  0.600511703
#11:       b   1  1.5117812           NA
#12:       b   2  0.3898432           NA
#13:       b   3 -0.6212406  0.950812202
#14:       b   4 -2.2146999  0.426794608
#15:       b   5         NA -0.815365744
#16:       b   6         NA -1.417970234
#17:       b   7         NA           NA
#18:       b   8  0.9438362           NA
#19:       b   9  0.8212212           NA
#20:       b  10  0.5939013  0.882528703

【讨论】：

不错的一个。您可能还想在代码中包含library(zoo)
嗯，是的，这是一个解决方案。然而，它只用 dplyrs group_by 替换 unlist(tapply()) ，用 dplyrs lag 函数替换我自己的 lag 函数。否则看起来一样。
@MrMax 是的，我同意。没有太大的不同，除了（也许）稍微更具可读性；我不清楚你的问题陈述，仍然很困惑。感谢您添加预期的输出，这似乎表明上述dplyr 代码不是正确的解决方案。
@MrMax 我已经根据set.seed(1) 的示例数据进行了更新。
不用担心@MrMax；我还添加了data.table 解决方案（应该比dplyr 方法更快）。

【解决方案3】：

这可能不是最优雅或可扩展的解决方案，但它确实提供了预期的结果：

df %>%
  group_by(subject) %>%
  mutate(n_values = 3 - is.na(lag(value, 1)) - is.na(lag(value, 2)) - is.na(lag(value, 3)),
         roll_mean = ifelse(
           n_values >= 2,
           (coalesce(lag(value), 0) + coalesce(lag(value, 2), 0) + coalesce(lag(value, 3), 0)) / n_values,
           NA)
  )

解释：这是一个dplyr 管道，它首先按主题分组，以便尊重组。接下来mutate中有两个计算值：

n_values 计算前 3 行中非 NA 值的个数，每个 NA 值等于 3 减 1。使用lag 访问前面的行。
roll_mean是有条件的，使用ifelse：如果n_values至少等于2，则可以计算均值。它将前 3 个值相加，使用 coalesce 将 NA 替换为 0。总和除以n_values 得到平均值。如果n_values < 2，则返回NA。

【讨论】：

这是一个有趣的解决方案。但是，作为中级 R 用户，我很难理解它是如何工作的。由于我正在寻找更简单的解决方案，可能不是预期的方向。
这可能有点令人费解，但由于我对rollapply 不是很熟悉，所以这是我最好的选择。我将编辑以在答案中添加一些对该过程的解释。