Arima预测模型（R语言）

转自 http://blog.csdn.net/desilting/article/details/39013825#comments

ARIMA(p，d，q)模型全称为差分自回归移动平均模型(Autoregressive Integrated Moving Average Model,简记ARIMA)，AR是自回归， p为自回归项； MA为移动平均，q为移动平均项数，d为时间序列成为平稳时所做的差分次数。

所谓ARIMA模型，是指将非平稳时间序列转化为平稳时间序列，然后将因变量仅对它的滞后值以及随机误差项的现值和滞后值进行回归所建立的模型。ARIMA模型根据原序列是否平稳以及回归中所含部分的不同，包括移动平均过程（MA）、自回归过程（AR）、自回归移动平均过程（ARMA）以及ARIMA过程。

下面我们选取了一些具有周期性（7天）的测试数据，通过ARIMA模型做一个简单的预测。

实现需要导入包

library(tseries)
library(forecast)
library(xts)
#install.packages("xts")安装包

source<- c(10930,10318,10595,10972,7706,6756,9092,10551,9722,10913,11151,8186,6422,
6337,11649,11652,10310,12043,7937,6476,9662,9570,9981,9331,9449,6773,6304,9355,10477,
10148,10395,11261,8713,7299,10424,10795,11069,11602,11427,9095,7707,10767,12136,12812,
12006,12528,10329,7818,11719,11683,12603,11495,13670,11337,10232,13261,13230,15535,
16837,19598,14823,11622,19391,18177,19994,14723,15694,13248,9543,12872,13101,15053,
12619,13749,10228,9725,14729,12518,14564,15085,14722,11999,9390,13481,14795,15845,
15271,14686,11054,10395,14775,14618,16029,15231,14246,12095,10473,15323,15381,14947)

测试数据的时间序列图，如下：

data<-xts(source,seq(as.POSIXct("2014-01-01"),len=length(data),by="day"))

plot(data)

通过以下代码：

data_diff1<-diff(data,differences=1)
plot(data_diff1) 
data_diff2<-diff(data,differences=2)
plot(data_diff2) 

可以看出一次差分后的时间序列在均值和方差上看起来像是平稳的，与二次差分的图形相差不大，随着时间推移，时间序列大致保持不变，因此设置差分项d=1。

接下来需要选择合适的ARIMA模型，即确定ARIMA(p,d,q)中合适的 p、q 值，我们通过R中的“acf()”和“pacf”函数来做判断。

1、查看自相关图

acf <- acf(data_diff1[!is.na(data_diff1)],lag.max=100,plot=FALSE)
plot(acf)

2、查看偏自相关图

 pacf <- pacf(data_diff1[!is.na(data_diff1)],lag.max=100,plot=FALSE) 
 plot(pacf)

ARMA(p，q)模型的ACF与PACF理论模式

模型	ACF	PACF
AR(p)	衰减趋于零（几何型或振荡型）	p阶后截尾
MA(q)	q阶后截尾	衰减趋于零（几何型或振荡型）
ARMA(p，q)	q阶后衰减趋于零（几何型或振荡型）	p阶后衰减趋于零（几何型或振荡型）

因此我们的arima模型为arima(7,1,0)

data.fit <- Arima(data,order=c(7,1,0), seasonal=list(order=c(1,1,0), period=7))
data.fit 

Series: data ARIMA(7,1,0)(1,1,0)[7] Coefficients:           ar1 ar2 ar3 ar4 ar5 ar6 ar7 sar1       -0.2829 -0.2128 -0.0180 0.0606 0.3164 0.0415 -0.0883 -0.5075 s.e. 0.1040 0.1076 0.1076 0.1048 0.1062 0.1066 0.1520 0.1449 sigma^2 estimated as 1589493: log likelihood=-789.34 AIC=1596.68 AICc=1598.88 BIC=1619.38

预测后一周的值：

forecast <- forecast.Arima(data.fit,h=7,level=c(99.5))
forecast
plot.forecast(forecast)

Forecast Lo 80 Hi 80 Lo 95 Hi 95 101 15206.61 13590.892 16822.33 12735.583 17677.64 102 14654.18 12665.935 16642.42 11613.422 17694.94 103 11589.60 9388.579 13790.62 8223.428 14955.77 104 10190.69 7742.040 12639.34 6445.803 13935.58 105 14362.50 11640.733 17084.26 10199.919 18525.07 106 14785.77 11620.577 17950.96 9945.025 19626.52 107 15501.02 12049.124 18952.91 10221.802 20780.23

tsdiag检验——

tsdiag(data.fit)

Acf检验说明：

      残差没有明显的自相关性

Ljung-Box测试显示：

      所有的P-value>0.5， 说明残差为白噪声