如何计算R中不规则降雨时间序列中的最大60分钟总数？答案

【问题标题】：How to calculate the maximum 60 minute total in an irregular rainfall timeseries in R?如何计算R中不规则降雨时间序列中的最大60分钟总数？
【发布时间】：2021-05-20 23:12:16
【问题描述】：

我有一个使用 0.5 毫米翻斗记录的持续不规则时间序列的降雨数据。

例如

Date and time	Rainfall
11/05/2021 11:05:17	0.5
11/05/2021 11:15:10	0.5
11/05/2021 11:20:04	0.5
11/05/2021 11:28:22	0.5
11/05/2021 11:33:25	0.5
11/05/2021 11:36:39	0.5
11/05/2021 11:39:50	0.5
11/05/2021 11:41:43	0.5
11/05/2021 11:43:35	0.5
11/05/2021 11:44:57	0.5
11/05/2021 11:47:02	0.5
11/05/2021 11:48:42	0.5
11/05/2021 11:53:04	0.5
11/05/2021 11:58:33	0.5
11/05/2021 12:01:27	0.5
11/05/2021 12:02:52	0.5
11/05/2021 12:07:35	0.5
11/05/2021 12:10:32	0.5
11/05/2021 12:12:55	0.5
11/05/2021 12:16:22	0.5
11/05/2021 12:17:45	0.5
11/05/2021 12:20:14	0.5
11/05/2021 12:22:26	0.5

在活动期间我希望能够计算：

总共滚动 60 分钟
活动总共最多 60 分钟（不是从上午 11 点到中午 12 点、中午 12 点到下午 1 点等的每小时总数）

即在上面的示例中：一个 60 分钟的时间段可以从 12:07:35 开始并包括回 11:15:10（60 分钟的时间段从小费的确切 hh:mm:ss 开始）。

了解最大总数有助于我们比较事件和预测。

到目前为止我的思考过程：

我一直在使用滞后函数计算水位上升率（每隔 5 分钟记录一次）。

River_Hour_RoR <- mutate(River, RoR = Stage - lag(Stage, n = 12))

我想也许我可以将 sum/cumulative-sum 函数与 lag 函数结合使用（类似于 WL 上升率），但我不确定如何指定 60 分钟的时间间隔。

任何关于如何做到这一点的想法或不同的方法将不胜感激！谢谢:)！

【问题讨论】：

标签： r lag cumulative-sum

【解决方案1】：

library(lubridate)
library(slider)


# convert to POSIXct datetime format
df1$Date.and.time = lubridate::mdy_hms(df1$Date.and.time)

# sum over prior 60 minutes using slider::slide_index
df1$hourly_total = slider::slide_index_dbl(df1$Rainfall, df1$Date.and.time, sum, .before = minutes(60))

df1[df1$hourly_total == max(df1$hourly_total),]
#         Date.and.time Rainfall hourly_total
#23 2021-11-05 12:22:26      0.5           10

或者，这里有一个dplyr 方法，我们将数据的累积总和加上一个滞后一小时的版本，其中降雨量为负。这将允许仅累积过去 60m 的总数。

library(dplyr)
bind_rows(
  df1,
  df1 %>% mutate(Date.and.time = Date.and.time + dhours(1),
                 Rainfall = -Rainfall)) %>%
  arrange(Date.and.time) %>%
  mutate(Rainfall_60m = cumsum(Rainfall))

如果我们将它输入 ggplot，我们可以直观地看到它是如何工作的：

... %>%
ggplot(aes(Date.and.time, Rainfall_60m)) + 
geom_step() +
geom_col(aes(y = Rainfall))

【讨论】：

【解决方案2】：

使用 data.table，仅对过去 60 分钟内的那些记录进行自联接：

library(data.table)
setDT(dat)

dat[, Datetime := as.POSIXct(Datetime , format="%d/%m/%Y %H:%M:%S", tz="UTC") ]
dat[, subthour := Datetime - as.difftime(1, units="hours") ]

dat[, Sum_Rainfall := 
  dat[
    dat,
    on = c("Datetime>=subthour", "Datetime<=Datetime"),
    sum(Rainfall),
    by=.EACHI
  ]$V1
]
dat[, subthour := NULL]

结果：

#               Datetime Rainfall Sum_Rainfall
# 1: 2021-05-11 11:05:17      0.5          0.5
# 2: 2021-05-11 11:15:10      0.5          1.0
# 3: 2021-05-11 11:20:04      0.5          1.5
# 4: 2021-05-11 11:28:22      0.5          2.0
# 5: 2021-05-11 11:33:25      0.5          2.5
# 6: 2021-05-11 11:36:39      0.5          3.0
# 7: 2021-05-11 11:39:50      0.5          3.5
# 8: 2021-05-11 11:41:43      0.5          4.0
# 9: 2021-05-11 11:43:35      0.5          4.5
#10: 2021-05-11 11:44:57      0.5          5.0
#11: 2021-05-11 11:47:02      0.5          5.5
#12: 2021-05-11 11:48:42      0.5          6.0
#13: 2021-05-11 11:53:04      0.5          6.5
#14: 2021-05-11 11:58:33      0.5          7.0
#15: 2021-05-11 12:01:27      0.5          7.5
#16: 2021-05-11 12:02:52      0.5          8.0
#17: 2021-05-11 12:07:35      0.5          8.0
#18: 2021-05-11 12:10:32      0.5          8.5
#19: 2021-05-11 12:12:55      0.5          9.0
#20: 2021-05-11 12:16:22      0.5          9.0
#21: 2021-05-11 12:17:45      0.5          9.5
#22: 2021-05-11 12:20:14      0.5          9.5
#23: 2021-05-11 12:22:26      0.5         10.0
#               Datetime Rainfall Sum_Rainfall

dat 在哪里：

dat <- structure(list(Datetime = c("11/05/2021 11:05:17", "11/05/2021 11:15:10", 
"11/05/2021 11:20:04", "11/05/2021 11:28:22", "11/05/2021 11:33:25", 
"11/05/2021 11:36:39", "11/05/2021 11:39:50", "11/05/2021 11:41:43", 
"11/05/2021 11:43:35", "11/05/2021 11:44:57", "11/05/2021 11:47:02", 
"11/05/2021 11:48:42", "11/05/2021 11:53:04", "11/05/2021 11:58:33", 
"11/05/2021 12:01:27", "11/05/2021 12:02:52", "11/05/2021 12:07:35", 
"11/05/2021 12:10:32", "11/05/2021 12:12:55", "11/05/2021 12:16:22", 
"11/05/2021 12:17:45", "11/05/2021 12:20:14", "11/05/2021 12:22:26"
), Rainfall = c(0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 
0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 
0.5)), class = "data.frame", row.names = c(NA, -23L))

【讨论】：

【解决方案3】：

将第一列转换为 POSIXct，然后使用findInterval 查找时间 60 分钟（3600 秒）后的行号。 w[i] 定义为数据帧第 i 行的总行数。使用允许宽度为向量 w 的运行 rollapplyr。 w 可能很有趣，因为它显示了每个总数中涉及的行数。然后我们可以使用 which.max 来获取总数最大的行。

library(zoo)

tt <- as.POSIXct(DF[[1]], format = "%m/%d/%Y %H:%M:%S")
w <- seq_along(tt) - findInterval(tt - 3600, tt)
DF2 <- transform(DF, total = rollapplyr(Rainfall, w, sum, fill = NA), w = w,
              check.names = FALSE)

head(DF2)
##         Date and time Rainfall total w
## 1 11/05/2021 11:05:17      0.5   0.5 1
## 2 11/05/2021 11:15:10      0.5   1.0 2
## 3 11/05/2021 11:20:04      0.5   1.5 3
#3 4 11/05/2021 11:28:22      0.5   2.0 4
## 5 11/05/2021 11:33:25      0.5   2.5 5
## 6 11/05/2021 11:36:39      0.5   3.0 6

DF2[which.max(DF2$total), ]
##          Date and time Rainfall total  w
## 23 11/05/2021 12:22:26      0.5    10 20

【讨论】：