时间序列模型在R软件中的实现

Posted 2021-05-02 深度学习与资源共享

本期主要是分享ARIMA模型在R语言中的实现，需要有一定的基础知识，如果有疑问请自行百度或者后台回复询问。

前期准备

✦R软件✔

✦ARIMA模型知识✔

R软件可以去官网自行下载，软件和程序包的安装过程中的问题请自行百度。

以下是安装的过程

软件下载

任选一个中国区的镜像

安装成功后任选一个中国区的镜像，然后加载程序包

左右滑动查看更多

请下载forecast，tseries，TSA这三个程序包，并使用library(forecast)加载程序包，后面类似。

（如果加载不了镜像可能是IE浏览器出现了问题，请重置IE浏览器再试一试,后面涉及到数据包中的函数，请一定要先加载数据包再进行后面的操作）

实现流程

实现的主要步骤：

实现过程

①将数据处理成时间序列

主要涉及的函数ts（）

t<-ts(c(数据),start=c(开始时间),frequency=1) 函数中的frequency的取值不同表示的时间间隔不同,frequency=12表示时间单位是年，每一个时间单位中有12个均匀间隔的观察值；frequency=4表明时间单位是年，每一个时间单位中有4个季节观察值；frequency=365表明时间单位是年，每一个时间单位中有365个日期观察值。

②判断时间序列是否平稳

主要涉及的函数ndiffs（）、diff（）、plot（）

首先，我们需要画图判断其平稳性。例如R中自带的数据Nile（尼罗河流量）,在命令窗口中输入Nile即可见。首先加载library(forecast)，library(tseries)画出Nile的图形plot(Nile)

ndiffs()函数可以帮助我们判断需要进行几次差分。

ndiffs(Nile)

dNile <- diff(Nile)#进行差分处理。上面显示需要进行一次差分                                         

plot(dNile)#查看一次差分后的图形

进行adf检验，如果结果显著，则认为序列是平稳的

adf.test(dNile)

看p-value 的值。P值要求小于给定的显著性水平,一般是0.05、0.01等,p越接近于0越好；一般认为p小于0.05就认为是显著的。

③ARMA模型的定阶（p，q值的选取）

通过ACF图和PACF图来选择备选模型。

Acf(dNile)#画出acf图

Pacf(dNile)#画出pacf图

在确定差分平稳后，需要判断p和q，这里定阶方法有很多，因为p和q的确定也很复杂，不是一下子就可以确定的。主要有这么几种(1)观察法，直接观察，如果acf在q+1阶突然截断，在q处截尾，则为ma(q)序列，同理，pacf在p处截尾则为ar（p）序列，否则为arma（p,q）序列，二者结合进一步判断(2)参数检验，利用数理统计检验高阶模型的新增加的参数是否近似为零，检验模型残差的相关特性等(3)信息准则，确定一个与模型阶数有关的准则，如AIC、BIC等，既考虑拟合效果接近程度，又考虑参数个数。实际中往往多种方法综合应用，选择最合适的p,d,q.

④模型拟合

涉及函数arima（），arima(ts,order=c(p,d,q))

fit <- arima(dNile, order=c(0,1,1)) #将拟合命名为fit

fit

accuracy(fit)#查看具体信息

MASE是相对最为重要的指标。当采用平均绝对误差（MAE）作为模型预测准确性的相对指标的话，如果MASE>1, 则表明对于样本之外的预测比基于样本本身的朴素预测要差。也即如果我们期望样本之外的资料能非常接近样本本身的话，当出现MASE>1的情况，则对样本之外的预测值是不可取的。

⑤模型显著性检验和评价

一般来说，一个模型如果合适，那模型的残差应该满足均值为零的正态分布，并且对于任意的滞后阶数，残差的自相关系数应该为零。模型的残差应该满足正态分布（即残差间没有关联）

qqnorm(fit$residuals)    #画出qq图判断是否服从正态分布

qqline(fit$residuals)

如果数据满足正态分布，则数据中的点会落在图中的线上。

结果现实较好的。

Box.test(fit$residuals, type="Ljung-Box")#检验残差的自相关系数是否都为零

查看是否通过显著性检验。结果显著。

⑥预测

直接用forecast(fit, 4)预测，后面的数字代表预测几年。并做图。

forecast(fit, 4)

plot(forecast(fit, 4), xlab="Year", ylab="Annual Flow")

蓝色是估计值，浅灰，深灰代表80%，95%的置信区间。

附

  完整版代码

library(forecast)

library(tseries)

#将生成图聚合到一起

par(mfrow=c(2,3))

plot(Nile)

#需要进行差分的次数

ndiffs(Nile)

dNile <- diff(Nile)                                              

plot(dNile)

adf.test(dNile)

#通过ACF图和PACF图来选择备选模型。

Acf(dNile)

Pacf(dNile)

#根据ACF图和PACF图确定参数p,q

fit <- arima(Nile, order=c(0,1,1))                                 

fit

accuracy(fit)

#通过残差检验模型是否合理

qqnorm(fit$residuals)     

qqline(fit$residuals)

Box.test(fit$residuals, type="Ljung-Box")

#  Forecasting with an ARIMA model

forecast(fit, 3)

plot(forecast(fit, 3), xlab="Year", ylab="Annual Flow")

看了那么多是不是觉得很复杂，其实是有函数可以进行自动预测的

ARIMA的自动预测

可以通过直接用forecast()包中的auto.arima()函数实现ARIMA模型的自动选取。例如：

library(forecast)

fit <- auto.arima(data)

fit

forecast(fit, 3)

accuracy(fit)

有兴趣的朋友可以自己动手操作一下哦~