数据分析课设（R）【实践-时间序列】

背景问题：

利用ARIMA乘积模型对航空旅客数据集进行分析一年的预测

数据集

直接导入R的数据集AirPassengers：这是1949-1960年每月旅客总数；数据集来源Box and Jenkins的文献。

#先将数据集导入

x=AirPassengers

#查看数据集：

print(x)

     Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

1949 112 118 132 129 121 135 148 148 136 119 104 118

1950 115 126 141 135 125 149 170 170 158 133 114 140

1951 145 150 178 163 172 178 199 199 184 162 146 166

1952 171 180 193 181 183 218 230 242 209 191 172 194

1953 196 196 236 235 229 243 264 272 237 211 180 201

1954 204 188 235 227 234 264 302 293 259 229 203 229

1955 242 233 267 269 270 315 364 347 312 274 237 278

1956 284 277 317 313 318 374 413 405 355 306 271 306

1957 315 301 356 348 355 422 465 467 404 347 305 336

1958 340 318 362 348 363 435 491 505 404 359 310 337

1959 360 342 406 396 420 472 548 559 463 407 362 405

1960 417 391 419 461 472 535 622 606 508 461 390 432

初步分析

#对x进行log处理

lx=log(x)

#一阶差分

dlx=diff(lx)

#滞后十二阶的差分

ddlx=diff(dlx,12)

#输出以上对应的图像

plot.ts(cbind(x,lx,dlx,ddlx),main="")

图片说明：

1.原数据

2.进行了对数处理 ->（降低尺度，稳定方差）

3.一阶差分的对数数据图

4.滞后12步差分的对数数据图 -> （每一年的月份存在一定的相关性）

分析：

首先原数据明显有均值增大的趋势，同时方差也在增大。该序列是不平稳的，然而平稳性是时间序列分析的前提条件，故我们需要对不平稳的序列进行处理将其转换成平稳的序列。

然后对对数进行一阶差分后，将原先的均值上升趋势消除，但是可以看到，方差略大，因此考虑12个月的年度因素，进行滞后12步差分，随后可以看见在1954-1958之间方差被有效地降低了。

基础分析验证

#设定打印是两行一列

par(mfrow=c(2,1))

#月度的一阶差分数据图

monthplot(dlx)

#月度的滞后十二阶差分数据图

monthplot(ddlx)

分析：

可以看到一阶的时候，均值不稳定，但是在滞后12步差分后，均值稳定到0上。也是验证了我们上述的说明。

#对滞后12步差分进行平稳性检验

adf.test(ddlx)

    Augmented Dickey-Fuller Test

data:  ddlx

Dickey-Fuller = -5.1993, Lag order = 5, p-value = 0.01

alternative hypothesis: stationary

Warning message:

In adf.test(ddlx) : p-value smaller than printed p-value

结果显示P<0.01因此差分后的数据可以认为是平稳的。

最后，滞后12步差分的数据看起来是比较平稳的，所以可以开始建模。

ACF与PACF

数据平稳后，需要对模型定阶，即确定p、q的阶数。

首先需要检查ACF和PACF。

#绘图和打印的时间序列的ACF和PACF

acf2(ddlx)

分析：

可以发现在每个季节成分上，ACF是在s=1上截尾，PACF是一个拖尾的样子，所以姑且考虑季节性成分是SMA(0,1,1)

反正不对的话，后面检验的时候也不会通过。

再来看非季节的成分：发现不管是ACF还是PACF好像是拖尾的样子。暂定它是ARMA(1,1)

尝试建模：乘积季节模型ARIMA

接下来要检验乘积季节模型ARIMA（1,1,1）X（0,1,1）

#建模ARIMA（1,1,1）X（0,1,1）

sarima(lx,1,1,1,0,1,1,12)

Coefficients:

         ar1      ma1     sma1

      0.1960  -0.5784  -0.5643

s.e.  0.2475   0.2132   0.0747

sigma^2 estimated as 0.001341:  log likelihood = 244.95,  aic = -481.9

结果发现AR才0.196，不显著，所以从季节因素中减少一个参数

比如ARIMA（1,1,0）X（0,1,1）

#建模ARIMA（1,1,0）X（0,1,1）

sarima(lx,1,1,0,0,1,1,12)

Coefficients:

          ar1     sma1

      -0.3395  -0.5619

s.e.   0.0822   0.0748

sigma^2 estimated as 0.001367:  log likelihood = 243.74,  aic = -481.49

我们发现果然AR绝对值提高了，但是会不会有更好的情况呢？

比如ARIMA（0,1,1）X（0,1,1）

#建模ARIMA（0,1,1）X（0,1,1）

sarima(lx,0,1,1,0,1,1,12)

Coefficients:

          ma1     sma1

      -0.4018  -0.5569

s.e.   0.0896   0.0731

sigma^2 estimated as 0.001348:  log likelihood = 244.7,  aic = -483.4

发现AIC变小了而最大似然变大了而且系数也更好了。

因此我们肯定选后一个模型。

于是我们得到拟合模型：

ARIMA（0,1,1）X（0,1,1）

1-B121-Bxt=1+ΘB12(1+θB)wt 

进一步对ARIMA（0,1,1）X（0,1,1）的plot残差情况进行分析：（也是sarima(lx,0,1,1,0,1,1,12)直接可得）

 

总体看起来还是比较平稳的，虽然也存在个别异常点的情况，比如图中左就有2.0附近超出边界的情况，以及最后一个图有个别点在拒绝的边缘（蓝色虚线）试探。

数据预测

最后进行一个数据预测：对未来一年的数据进行月份预测。

#对sarima(lx,0,1,1,0,1,1,12)进一步预测后面12个月的情况

sarima.for(lx,12,0,1,1,0,1,1,12)

$pred

          Jan      Feb      Mar      Apr      May      Jun      Jul      Aug      Sep      Oct      Nov      Dec

1961 6.110186 6.053775 6.171715 6.199300 6.232556 6.368779 6.507294 6.502906 6.324698 6.209008 6.063487 6.168025

$se

            Jan        Feb        Mar        Apr        May        Jun        Jul        Aug        Sep        Oct

1961 0.03671562 0.04278291 0.04809072 0.05286830 0.05724856 0.06131670 0.06513124 0.06873441 0.07215787 0.07542612

            Nov        Dec

1961 0.07855851 0.08157070

 

其中阴影范围就是可能震荡的区间。