在R中基于日期条件筛选：匹配降水同位素与降水事件数据

核心思路 #

为每个同位素样本划定对应的降水观测窗口，统计窗口内的降水事件数量，仅保留对应单次降水事件的样本，最终关联降水日期、降水量和同位素值。

实现代码 #

library(dplyr)

# 原始数据（用户提供）
start <- as.Date('2017/01/01')
len <- 21
dates <- seq(start, by = "day", length.out = len)
prcp <- c(0, 1.0, 2.0, 0, 1.0, 0, 0,
          0, 0, 0, 0, 0, 1.75, 2.0,
          0, 0, 0, 0, 0, 0, 0)
samp <- c(NA, NA, NA, NA, -15.0, NA, NA,
          NA, NA, NA, NA, NA, -12.0, NA,
          NA, NA, NA, NA, NA, -20, NA) 
raw <- data.frame(dates, prcp, samp)

# 1. 提取所有有效样本（非NA的samp）
sample_points <- raw %>%
  filter(!is.na(samp)) %>%
  select(sample_date = dates, isotopes = samp) %>%
  mutate(prev_sample_date = lag(sample_date, default = start - 1)) # 第一个样本的前序日期设为数据起始前一天

# 2. 遍历每个样本，匹配对应的降水事件
result_list <- list()
for(i in 1:nrow(sample_points)){
  # 确定当前样本的降水观测窗口：上一个样本次日 到 当前样本日期
  window_start <- sample_points$prev_sample_date[i] + 1
  window_end <- sample_points$sample_date[i]
  
  # 提取窗口内的降水数据
  window_prcp <- raw %>%
    filter(dates >= window_start & dates <= window_end) %>%
    filter(prcp > 0)
  
  # 仅保留窗口内只有单次降水的情况
  if(nrow(window_prcp) == 1){
    result_list[[i]] <- data.frame(
      dates_out = window_prcp$dates,
      prcp_out = window_prcp$prcp,
      samp_out = sample_points$isotopes[i]
    )
  }
}

# 3. 合并结果并整理
final <- bind_rows(result_list)
print(final)

代码解释 #

• 提取样本点：筛选出所有有同位素测量值的日期，同时为每个样本标记上一个样本的日期，用于划分专属降水观测窗口。
• 划分降水窗口：每个样本的观测窗口设定为上一个样本采集日的次日到当前样本采集日，确保窗口对应唯一的样本周期。
• 筛选有效样本：统计窗口内非零降水的天数，仅保留天数为1的情况，将降水事件与对应同位素值关联。
• 合并结果：把符合条件的匹配结果整合成最终数据框。

运行结果 #

执行代码后会得到与期望一致的结果：

dates_out prcp_out samp_out
1 2017-01-13     1.75    -12.0
2 2017-01-14     2.00    -20.0

替代方案建议 #

如果数据量较大，可使用purrr包替代循环，实现更高效的批量处理：

library(purrr)

final <- sample_points %>%
  mutate(
    window_data = map2(prev_sample_date + 1, sample_date, function(start, end){
      raw %>% filter(dates >= start & dates <= end, prcp > 0)
    }),
    valid = map_lgl(window_data, ~nrow(.) == 1)
  ) %>%
  filter(valid) %>%
  transmute(
    dates_out = map_chr(window_data, ~.$dates),
    prcp_out = map_dbl(window_data, ~.$prcp),
    samp_out = isotopes
  )

内容的提问来源于stack exchange，提问作者Philip Moffatt