R_非季节性ARIMA模型

install.package(‘forecast’) library(readx) library(forecast) data <- read_xlsx(‘a.xlsx’) ts <- ts(data$sales,frequency=12,start = c(2012,8))

$$ \[\begin{equation} \begin{split} &自回归模型：序列中的每一个值都可以由它之前的p个值的线性组合来表示，具体表示为 \\ &AR(p):Y_t = c + \beta_1Y_{t-1}+ \beta_2Y_{t-2}+ ...+ \beta_pY_{t-p} +\varepsilon_t \\ &移动平均模型：时序中的每个值都可以用之前的q个残差的线性组合来表示，具体表示为：\\ &MA(q):Y_t = c + \varepsilon_{t} + \theta_1\varepsilon_{t-1}+ \theta_2\varepsilon_{t-2}+...+ \theta_q\varepsilon_{t-q}\\ &ARMA(p,q):Y_t = \beta_0 + \beta_1Y_{t-1}+ \beta_2Y_{t-2}+ ...+ \beta_1Y_{t-p} +\varepsilon_t -\theta_1\varepsilon_{t-1}+ \theta_1\varepsilon_{t-1}+...+ \theta_1\varepsilon_{t-q} \\ &\\ \end{split} \end{equation}\] $$

nhtemp

Time Series:
Start = 1912 
End = 1971 
Frequency = 1 
 [1] 49.9 52.3 49.4 51.1 49.4 47.9 49.8 50.9 49.3 51.9 50.8 49.6 49.3 50.6 48.4 50.7 50.9 50.6
[19] 51.5 52.8 51.8 51.1 49.8 50.2 50.4 51.6 51.8 50.9 48.8 51.7 51.0 50.6 51.7 51.5 52.1 51.3
[37] 51.0 54.0 51.4 52.7 53.1 54.6 52.0 52.0 50.9 52.6 50.2 52.6 51.6 51.9 50.5 50.9 51.7 51.4
[55] 51.7 50.8 51.9 51.8 51.9 53.0

平稳性检验(单位根) 和 差分变换

$$ \[\begin{equation} \begin{split} &平稳就是指序列的统计性质并不会随着时间的推移而改变，例如Y_t的均值和方差都是恒定的。\\ &另外，对任意滞后阶数k，此序列的自相关性不发生改变。\\ &时间序列的平稳性一般可以通过图形来直接判断，观察是否有明显的趋势向和不规则波动项\\ &也可使用单位根（ADF）统计检验来验证平稳性，原假设为该序列存在单位根，如果经过原假设不能被拒绝的话，那么我们就可以认为该序列为非平稳时间序列。\\ &R中的ADF检验可使用tseries包中的adf.test()函数来实现： \end{split} \end{equation}\] $$

adf.test(nhtemp)

    Augmented Dickey-Fuller Test

data:  nhtemp
Dickey-Fuller = -3.2773, Lag order = 3, p-value = 0.08376
alternative hypothesis: stationary

# 备择假设为'stationary',即平稳时序,零假设也就自然为非平稳时序
# 已知P值为0.08>0.05,无法在95%的置信度下拒绝原假设,所以该序列不满足平稳性假设

ndiffs(nhtemp)

[1] 1

dnhtemp <- diff(nhtemp,1)
adf.test(dnhtemp)

p-value smaller than printed p-value

    Augmented Dickey-Fuller Test

data:  dnhtemp
Dickey-Fuller = -4.6366, Lag order = 3, p-value = 0.01
alternative hypothesis: stationary

# 根据P值可拒绝原假设,一次差分后实现了平稳性

par(mfrow=c(1,2))
plot(nhtemp)
plot(diff(nhtemp,1))

模型定阶

par(mfrow=c(1,2))
acf(dnhtemp)
Pacf(dnhtemp)

以	此准则为指导观察ACF与PACF图，可以发现，dnhtemp数据集在滞后项为一阶时有非常明显的自相关，而当滞后阶数逐渐增加时，偏相关基本呈现逐渐减小的趋势（至少没有超过置信区间），所以可以考虑arima(0,1,1)模型。
|	模型 | ACF | PACF |
|A	RIMA(p,d,0)|逐渐减小到0 |在P阶后减小到0 |
|A	RIMA(0,d,q)|q阶后减小到0|逐渐减小到0 |
|A	RIMA(p,d,q)|逐渐减小到0 |逐渐减小到0 |

#forecast包中的auto.arima()函数可以帮助我们自行选取最优参数值，
# 这比观察ACF图自行确定要便捷的多。
model_fit <- auto.arima(nhtemp, seasonal = FALSE)
model_fit

Series: nhtemp 
ARIMA(0,1,1) 

Coefficients:
          ma1
      -0.7983
s.e.   0.0956

sigma^2 estimated as 1.313:  log likelihood=-91.76
AIC=187.52   AICc=187.73   BIC=191.67

#模型验证 模型的残差应该满足独立正态分布

model_fit <- arima(nhtemp,order=c(0,1,1))
#accuracy(model_fit)
summary(model_fit)

Call:
arima(x = nhtemp, order = c(0, 1, 1))

Coefficients:
          ma1
      -0.7983
s.e.   0.0956

sigma^2 estimated as 1.291:  log likelihood = -91.76,  aic = 187.52

Training set error measures:
                    ME     RMSE       MAE       MPE     MAPE      MASE         ACF1
Training set 0.1285716 1.126722 0.8952077 0.2080785 1.749865 0.7513123 -0.008258617

残差正态性

qqnorm(model_fit$residuals)
qqline(model_fit$residuals)

残差间相关性

acf(model_fit$residuals)

pacf(model_fit$residuals)

Box.test(model_fit$residuals,type='Ljung-Box')

    Box-Ljung test

data:  model_fit$residuals
X-squared = 0.0043004, df = 1, p-value = 0.9477

# P值为0.95,说明模型的残差没有通过显著性检验(即可以认为残差的自相关系数为零),ARIMA模型能较好地拟合本数据

#模型预测

# 以下值分别为预测值及80%,95%的置信区间值
forecast(model_fit,5)

plot(forecast(model_fit,5))