🍉CSDN小墨&晓末:https://blog.csdn.net/jd1813346972
个人介绍: 研一｜统计学｜干货分享
擅长Python、Matlab、R等主流编程软件
累计十余项国家级比赛奖项，参与研究经费10w、40w级横向
文章目录
- 1 目的
- 2 原序列差分处理
- 3 差分后序列平稳性检验
- 4 差分后序列白噪声检验
- 5 ARIMA乘法模型建立
- 6 ARIMA乘法模型拟合
- 7 ARIMA乘法模型显著性检验
- 8 ARIMA乘法疏系数模型
- 9 ARIMA乘法疏系数模型显著性检验
- 10 ARIMA乘法疏稀疏模型预测
1 目的
该篇文章主要展示针对时序进行ARIMA乘法模型建模，并根据实际情况进行改进。案例数据同 时间序列分析实战（三）：时序因素分解法：某欧洲小镇1963年1月至1976年12月每月旅馆入住的房间数构成时间序列 x t x_t xt，我们在 时间序列分析实战（五）：ARIMA加法（疏系数）模型建模利用此数据构建了ARIMA加法模型。现考虑演示根据该数据建立ARIMA乘法模型。
2 原序列差分处理
从 时间序列分析实战（三）：时序因素分解法一文中可知，该序列具有趋势和季节效应，进行1阶差分提取趋势效应，12步差分提取季节效应。
运行程序：
#对原数据进行1阶12步差分
y=diff(diff(data1,12))
plot(y,sub='图1 入住房间数差分后序列时序图')
运行结果：
3 差分后序列平稳性检验
运行程序：
#差分后序列平稳性检验
library(aTSA)
adf.test(y)
运行结果：
## Augmented Dickey-Fuller Test
## alternative: stationary
##
## Type 1: no drift no trend
## lag ADF p.value
## [1,] 0 -19.56 0.01
## [2,] 1 -11.01 0.01
## [3,] 2 -10.63 0.01
## [4,] 3 -9.08 0.01
## [5,] 4 -10.57 0.01
## Type 2: with drift no trend
## lag ADF p.value
## [1,] 0 -19.50 0.01
## [2,] 1 -10.98 0.01
## [3,] 2 -10.60 0.01
## [4,] 3 -9.05 0.01
## [5,] 4 -10.53 0.01
## Type 3: with drift and trend
## lag ADF p.value
## [1,] 0 -19.44 0.01
## [2,] 1 -10.94 0.01
## [3,] 2 -10.56 0.01
## [4,] 3 -9.01 0.01
## [5,] 4 -10.49 0.01
## ----
## Note: in fact, p.value = 0.01 means p.value <= 0.01
4 差分后序列白噪声检验
运行程序：
#差分后序列纯随机性检验
for(k in 1:2) print(Box.test(y,lag=6*k,type="Ljung-Box"))
运行结果：
##
## Box-Ljung test
##
## data: y
## X-squared = 56.87, df = 6, p-value = 1.941e-10
##
##
## Box-Ljung test
##
## data: y
## X-squared = 76.751, df = 12, p-value = 1.713e-11
通过1阶12步差分后序列时序图(图1)显示差分后的序列没有明显趋势和周期特征了，ADF检验结果显示差分后序列平稳，纯随机性检验结果显示差分后序列为非白噪声序列，适合建模。
5 ARIMA乘法模型建立
运行程序：
par(mfrow=c(2,1))
acf(y,lag=36)
title(sub="图2 入住房间数差分后序列自相关图(ACF)")
pacf(y,lag=36)
title(sub="图3 入住房间数差分后序列偏自相关图(PACF)")
运行结果：
为考虑季节相关特征，这时考察延迟12阶、24阶等以周期长度为单位的自相关系数和偏自相关系数的特征。图2和图3显示,自相关图表现拖尾性质，偏自相关图表现截尾性质，若把偏自相关系数看作1周期截尾，结合前面的差分信息，经过多次尝试后，尝试拟合ARIMA乘法模型:： A R I M A ( 3 , 1 , 2 ) × ( 1 , 1 , 0 ) 12 ARIMA(3,1,2)×(1,1,0)_{12} ARIMA(3,1,2)×(1,1,0)12
6 ARIMA乘法模型拟合
运行程序：
fit4=arima(data1,order = c(3,1,2),seasonal = list(order=c(1,1,0),period=12))
fit4
运行结果：
##
## Call:
## arima(x = data1, order = c(3, 1, 2), seasonal = list(order = c(1, 1, 0), period = 12))
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 sar1
## 0.8229 0.0030 -0.3107 -1.7905 0.7905 -0.4073
## s.e. 0.1155 0.1036 0.0851 0.1110 0.1102 0.0816
##
## sigma^2 estimated as 213.2: log likelihood = -640.63, aic = 1295.25
7 ARIMA乘法模型显著性检验
运行程序：
#模型显著性检验
ts.diag(fit4)
title(sub="图4 ARIMA(3,1,2)×(1,1,0)_12模型显著性检验")
运行结果：
此时模型 ϕ 2 \phi _2 ϕ2估计值在两倍标准差内，考虑建立ARIMA乘法疏系数模型： A R I M A ( ( 1 , 3 ) , 1 , 2 ) × ( 1 , 1 , 0 ) 12 ARIMA((1,3),1,2)×(1,1,0)_{12} ARIMA((1,3),1,2)×(1,1,0)12
8 ARIMA乘法疏系数模型
运行程序：
#拟合乘法ARIMA模型
fit5=arima(data1,order = c(3,1,2),seasonal = list(order=c(1,1,0),period=12),transform.pars = F,fixed=c(NA,0,NA,NA,NA,NA))
fit5
运行结果：
##
## Call:
## arima(x = data1, order = c(3, 1, 2), seasonal = list(order = c(1, 1, 0), period = 12),
## transform.pars = F, fixed = c(NA, 0, NA, NA, NA, NA))
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 sar1
## 0.8263 0 -0.3090 -1.7930 0.793 -0.4072
## s.e. 0.0941 0 0.0611 0.1078 0.107 0.0807
##
## sigma^2 estimated as 213.2: log likelihood = -640.63, aic = 1293.25
9 ARIMA乘法疏系数模型显著性检验
运行程序：
#模型显著性检验
ts.diag(fit5)
title(sub="图5 ARIMA((1,3),1,2)模型显著性检验")
运行结果：
此时，模型参数估计值在2倍标准差外，且残差为白噪声序列，aic=1293.25。最终 A R I M A ( ( 1 , 3 ) , 1 , 2 ) × ( 1 , 1 , 0 ) 12 ARIMA((1,3),1,2)×(1,1,0)_{12} ARIMA((1,3),1,2)×(1,1,0)12模型为：
∇ ∇ 12 x t = 1 − 0.79 B − 0.79 B 2 1 − 0.83 B + 0.31 B 3 ( 1 − 0.41 B 1 2 ) ϵ t , V a r ( ϵ t ) = 213.2 \nabla\nabla_{12}x_t=\frac{1-0.79B-0.79B^2}{1-0.83B+0.31B^3}(1-0.41B^12)\epsilon_t,Var(\epsilon_t)=213.2 ∇∇12xt=1−0.83B+0.31B31−0.79B−0.79B2(1−0.41B12)ϵt,Var(ϵt)=213.2
综上所述，从AIC的角度来看，疏系数的ARIMA乘法模型的AIC值更小，为1293.25，因此在该序列的拟合模型中，选用ARIMA乘法模型 A R I M A ( ( 1 , 3 ) , 1 , 2 ) × ( 1 , 1 , 0 ) 12 ARIMA((1,3),1,2)×(1,1,0)_{12} ARIMA((1,3),1,2)×(1,1,0)12模型来进行预测。
10 ARIMA乘法疏稀疏模型预测
运行程序：
library(forecast)
fore2=forecast::forecast(fit5,h=36)
fore2$mean
运行结果：
## Jan Feb Mar Apr May Jun Jul
## 1977 827.6433 768.2262 771.2960 860.6892 838.8292 963.9901 1128.1026
## 1978 853.5056 785.6367 792.7410 886.3393 869.2913 985.3445 1154.0686
## 1979 875.9526 811.5455 817.0429 908.9636 889.9786 1009.7484 1176.5900
## Aug Sep Oct Nov Dec
## 1977 1159.3063 892.6915 897.5451 778.5312 884.5493
## 1978 1178.2927 915.9224 915.4218 804.8587 914.4704
## 1979 1203.6458 939.5356 941.2077 827.1998 935.3489
运行程序：
plot(fore2,lty=2,sub='图5 入住房间数序列疏系数的ARIMA乘法模型预测效果图')
lines(fore2$fitted,col=4)
运行结果：
根据结果显示，疏稀疏的ARIMA乘法模型有较好拟合效果。
