如何通过Quadratic/parabolic interpolation插补加州新冠病例率异常数据？

你观测到的圣诞节前后的低谷是典型的假期检测上报延迟伪影，通过二次插值修正的具体实现如下：

实现步骤 #

• 第一步：定位伪影区间，排除2020-12-12至2021-01-13的异常数据，提取区间前后未受假期影响的可靠数据点用于拟合二次抛物线
• 第二步：将日期格式转换为数值型天数（以序列起始日为基准的偏移天数），满足多项式拟合的输入要求
• 第三步：用可靠数据点拟合二阶多项式模型
• 第四步：用拟合好的模型预测伪影区间内所有日期的病例数值，替换原始数据中对应位置的异常值
• 第五步：绘制修正前后的折线图，验证插值效果

完整可运行R代码 #

# 加载所需包
library(ggplot2)

# 原始数据
x <- seq.Date(as.Date("2020-09-01"), as.Date("2021-03-01"), by=1)
y <- c(9.36,9.16,9.05,8.88,8.76,8.65,7.94,7.81,7.65,7.5,7.47,7.5,7.52,8.19,
       8.03,8.1,8.12,8.14,8.19,8.24,8.21,8.19,8.19,8.22,8.24,8.2,8.16,8.14,
       8.16,8.14,8.25,8.19,8.3,8.36,8.45,8.43,8.42,8.44,8.51,8.63,8.62,8.63,
       8.66,8.69,8.73,8.81,8.79,8.9,9.15,9.46,9.67,9.78,10.07,10.19,10.32,10.48,
       10.52,10.69,10.93,11.27,11.68,12.4,13.45,14.66,15.92,17.09,18.15,18.85,
       19.04,19.98,20.93,21.69,22.89,24.28,25.52,26.78,29.08,31.29,33.62,35.34,
       37.11,37.95,38.35,39.59,40.82,42.06,39.44,39.63,41.69,43.73,47.78,52.64,
       57.16,65.24,70.15,72.29,73.01,76.01,78.53,81.46,84.64,87.58,89.86,90.79,
       93.81,96.47,98.05,99.48,100.07,99.73,99.65,99.36,99.52,99.32,92.84,82.53,
       84.34,86.33,89.6,92.99,96.15,99.42,101.56,102.45,102.72,103.63,102.26,101,
       104.48,112.58,109.57,106.79,100.29,94.47,88.73,83.79,79.33,76.19,74.25,
       67.69,63.86,60.59,57.27,54.07,51.8,50.69,49.49,46.01,42.54,39.79,37.28,
       36.11,35.29,33.53,31.66,30.16,28.58,27.37,26.13,25.13,23.06,21.33,19.92,
       18.65,17.51,16.71,16.13,14.63,13.89,13.03,12.27,11.56,11.1,10.79,10.63,
       10.07,9.63,9.28,8.98,8.77,8.61,8.25)
df <- data.frame(date = x, raw_y = y)

# 定义伪影区间，可根据实际情况调整
artifact_start <- as.Date("2020-12-12")
artifact_end <- as.Date("2021-01-13")

# 提取拟合用的可靠数据（排除伪影区间）
fit_df <- df[df$date < artifact_start | df$date > artifact_end, ]
# 转换日期为数值型（距序列起始日的天数）
fit_df$day_num <- as.numeric(fit_df$date - min(fit_df$date))
df$day_num <- as.numeric(df$date - min(df$date))

# 拟合二次多项式模型
quad_model <- lm(raw_y ~ poly(day_num, 2), data = fit_df)

# 预测伪影区间的修正值
artifact_idx <- which(df$date >= artifact_start & df$date <= artifact_end)
df$corrected_y <- df$raw_y
df$corrected_y[artifact_idx] <- predict(quad_model, newdata = df[artifact_idx, ])

# 绘制修正前后对比图
ggplot(df) +
  geom_line(aes(x = date, y = raw_y, color = "原始数据"), linewidth = 0.8) +
  geom_line(aes(x = date, y = corrected_y, color = "修正后数据"), linewidth = 0.8) +
  labs(x = "日期", y = "每10万居民新冠病例数", color = "数据类型") +
  theme_minimal()

结果说明 #

修正后的数据会消除圣诞节前后的异常低谷，匹配疫情整体上升-峰值-下降的二次趋势，更接近真实的病例增长情况。如果需要调整拟合范围，仅修改artifact_start和artifact_end的取值即可。

内容的提问来源于stack exchange，提问作者David Crow