如何使用R forecast包纳入气温水温变量预测鱼类数量并可视化

1. 前置说明 #

要将气温、水温作为外部特征纳入时间序列模型，我们使用带外部回归项的ARIMA模型（ARIMAX），对应forecast包中auto.arima的xreg参数即可实现。需要注意：预测未来值时，也需要提供对应预测周期的外部变量（气温、水温）数值，若没有实际监测的未来值，可先对两个温度变量单独做时间序列预测得到未来值。

另外你的数据存在两处需要预处理的问题：

• 包含221、222两个区域的观测，建议分开建模效果更优
• airtemp_f列第39行数值为174.5，属于明显异常值，我们先使用前后值的均值替换

2. 完整实现代码 #

# 加载所需包
library(forecast)
library(ggplot2)
library(xts)

# 预处理异常值
fish$airtemp_f[39] <- mean(c(fish$airtemp_f[38], fish$airtemp_f[40]))

# 示例：拿221区域的数据演示，222区域操作逻辑完全一致
fish_221 <- fish[fish$district == "221", ]
# 按日期排序保证时间序列顺序
fish_221 <- fish_221[order(fish_221$date), ]

# 构造历史时间序列和外部回归矩阵
count_ts <- ts(fish_221$fishcount, frequency = 10) # 每年10个周期，设置季节频率
xreg_history <- model.matrix(~ airtemp_f + watertemp_f, data = fish_221)[, -1] # 去掉截距项

# 训练ARIMAX模型
arimax_model <- auto.arima(count_ts, xreg = xreg_history)

# 构造未来3期的外部变量值：先对两个温度变量做预测
# 气温预测
air_ts <- ts(fish_221$airtemp_f, frequency = 10)
air_forecast <- forecast(auto.arima(air_ts), h = 3)
# 水温预测
water_ts <- ts(fish_221$watertemp_f, frequency = 10)
water_forecast <- forecast(auto.arima(water_ts), h = 3)
# 拼接成预测用的xreg矩阵
xreg_future <- cbind(airtemp_f = air_forecast$mean, watertemp_f = water_forecast$mean)

# 执行鱼量预测，h=3对应未来3个周期
fish_forecast <- forecast(arimax_model, xreg = xreg_future, level = 0.95, h = 3)
# 查看预测结果
print(fish_forecast)

# 用ggplot绘制预测结果
# 方法1：直接用forecast包自带的autoplot方法（基于ggplot2）
autoplot(fish_forecast) +
  labs(title = "221区域未来3周期鱼类数量预测结果",
       x = "时间周期",
       y = "鱼类数量",
       caption = "阴影部分为95%置信区间") +
  theme_minimal()

# 方法2：自定义ggplot绘制（适合需要灵活调整样式的场景）
# 先把历史数据和预测结果整理成数据框
history_df <- data.frame(
  period = 1:length(count_ts),
  value = as.numeric(count_ts),
  type = "历史观测值"
)
forecast_df <- data.frame(
  period = (length(count_ts)+1):(length(count_ts)+3),
  value = as.numeric(fish_forecast$mean),
  lower = as.numeric(fish_forecast$lower),
  upper = as.numeric(fish_forecast$upper),
  type = "预测值"
)
ggplot() +
  geom_line(data = history_df, aes(x = period, y = value, color = type)) +
  geom_line(data = forecast_df, aes(x = period, y = value, color = type)) +
  geom_ribbon(data = forecast_df, aes(x = period, ymin = lower, ymax = upper), 
              alpha = 0.2, fill = "blue") +
  scale_color_manual(values = c("历史观测值" = "black", "预测值" = "blue")) +
  labs(title = "221区域未来3周期鱼类数量预测结果",
       x = "观测周期",
       y = "鱼类数量",
       color = "数据类型",
       caption = "阴影部分为95%置信区间") +
  theme_minimal()

3. 关键说明 #

• 如果有实际监测的未来气温、水温数值，直接替换xreg_future的数值即可，不需要单独预测温度变量，预测精度会更高
• 222区域的预测只需要把筛选条件改成fish$district == "222"，其余代码逻辑完全一致
• 如果不需要按区域拆分建模，去掉区域筛选步骤直接用全量数据建模即可

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