基于WIN10的64位系统演示
一、写在前面
这一期，我们使用R进行SARIMA模型的构建。
同样，这里使用这个数据：
《PLoS One》2015年一篇题目为《Comparison of Two Hybrid Models for Forecasting the Incidence of Hemorrhagic Fever with Renal Syndrome in Jiangsu Province, China》文章的公开数据做演示。数据为江苏省2004年1月至2012年12月肾综合症出血热月发病率。运用2004年1月至2011年12月的数据预测2012年12个月的发病率数据。
R语言，这个大家都很熟悉吧。说实话，我不熟悉。
这里用的版本是：R-4.3.1，搭配RStudio食用。
二、R建立SARIMA实战
（1）导入数据
打开Rstudio，如下图操作，导入数据：
（2）单位根（ADF）检验
library(tseries)
# 假设你的时间序列数据存储在变量data$incidence中
test_result <- adf.test(data$incidence, alternative = "stationary")
# 打印测试结果
print(test_result)结果显示是平稳的，6666。那就平稳吧：
（3）设置为时间序列格式
# 将时间列转换为日期格式
time_series <- ts(data$incidence, frequency=12, start=c(2004, 1))（4）差分
# 进行一次差分
first_difference <- diff(time_series, differences = 1)
plot(first_difference, main="一次差分后的时间序列")
# 进行季节性差分
seasonal_difference <- diff(time_series, lag = 12)
plot(seasonal_difference, main="季节性差分后的时间序列")如图：
看起来，一次差分也得了。
（5）自相关图和偏相关图
绘制自相关图 (ACF)
acf(time_series, main="自相关图 (ACF)")
# 绘制偏自相关图 (PACF)
pacf(time_series, main="偏自相关图 (PACF)")如图，有点丑：
（6）建模
（6.1）数据拆分
# 划分训练集和验证集
train_series <- window(time_series, start=c(2004,1), end=c(2011,12))
validation_series <- window(time_series, start=c(2012,1), end=c(2012,12))（6.2）搭建SARIMA
看代码，自行体会：
案例：SARIMA(0,1,1)(0,1,1)12：
sarima_model <- Arima(train_series, order=c(0,1,1), seasonal=list(order=c(0,1,1), period=12))解读：这里，order=c(0,1,1)定义了非季节部分的阶数，而seasonal=list(order=c(0,1,1), period=12)定义了季节部分的阶数和季节周期（在这种情况下为12）。
（6.3）看模型参数
# 显示模型摘要
summary(sarima_model)如图：
解读如下：
参数没有统计学差异，这模型不得。大家自己试了，我继续。
（6）预测
（6.1）拟合数据
# 获取拟合数据
fitted_values <- fitted(sarima_model)
# 保存拟合数据到CSV文件
write.csv(data.frame(time=as.character(time(fitted_values)), fitted_values=fitted_values), file="fitted_values.csv", row.names=FALSE)
# 打印消息
print("拟合数据已保存到fitted_values.csv文件中")知道存在哪里不，看R的工作路径：
# 获取当前工作目录
current_working_directory <- getwd()
# 打印当前工作目录
print(current_working_directory)
（6.1）预测数据
# 预测未来12个月的数据
future_forecast <- forecast(sarima_model, h=12)
# 将预测结果与时间戳合并为数据框
forecast_data <- data.frame(
time = time(future_forecast$mean),
forecast_values = as.numeric(future_forecast$mean)
)
# 保存预测结果到CSV文件
write.csv(forecast_data, file="forecast_values.csv", row.names=FALSE)
# 打印消息
print("预测数据已保存到forecast_values.csv文件中")收工！！
四、数据
链接：https://pan.baidu.com/s/1qOpPi9pfzKR8TVmpOZaZcg?pwd=tc2z
提取码：tc2z
