1.ARMA过程
,
, 其中
为白噪声，且
。
2.过程平稳性、可逆性的判定
过程是平稳的当且仅当方程
解的模大于1.
过程是可逆的当且仅当方程
解的模大于1.
3.模型阶的识别
若时间序列是平稳的，则有三种方法来确定阶数的问题。
(1)第一种方法: acf图和pacf图
--------------------------------------
AR(p) MA(q) ARMA(p,q)
acf 图 拖尾 q阶截尾前q个无规律，其后拖尾
pacf图 p阶截尾 拖尾前p个无规律，其后拖尾
------------------------------------------
该方法为识别AR模型和MA模型提供了有效的工具，但对于ARMA模型来说，仍然非常困难。
(2)第二种方法: eacf图
Tsay和Tiao(1984年)提出了利用推广的自相关系数(eacf)来确定ARMA过程的阶数。eacf的推导比较复杂。这里根据R让软件说明eacf表如何确定ARMA过程的阶数。该表的行对应的是AR的阶，列对应的是MA的阶。该表的主要特征是：包含由“0”组成的三角形，并且这个三角形左上角顶点位于阶(p,q)处，使用这样的特征来识别ARMA过程的阶。对白噪声求eacf表
library(TSA)
x
eacf(x)＃给出白噪声的eacf表格
运行结果
－－－－－－－－－－－－－－
AR/MA
0 1 2 3 4 5 6 7 8 9 10 11 12 13
0 o o o o o o o o o o o o o o
1 x o o o o o o o o o o o o o
2 x x o o o o o o o o o o o o
3 x x x o o o o o o o o o o o
4 x x x o o o o o o o o o o o
5 x x x x x o o o o o o o o o
6 x x x x x x o o o o o o o o
7 x x x x o o x o o o o o o o
－－－－－－－－－－－－－－对AR(2)过程求eacf表格
x
eacf(x)
运行结果
－－－－－－－－－－－－－－
AR/MA
0 1 2 3 4 5 6 7 8 9 10 11 12 13
0 x x x x x x x x x x x x x x
1 x x o o o o o o o o x x o o
2 o o o o o o o o o o o x x o
3 x o o o o o o o o o o x o o
4 x x o o o o o o o o o o o o
5 x o x o x x o o o o o o o o
6 x o x x x x o o o o o o o o
7 x x o x x x o o o o o o o o
－－－－－－－－－－－－－－
对MA(2)过程求eacf表格
x
eacf(x)
运行结果
－－－－－－－－－－－－－－
AR/MA
0 1 2 3 4 5 6 7 8 9 10 11 12 13
0 x x o o o o o o o o o o o o
1 x x x o o o o o o o o o o o
2 x x x x o o o o o o o o o o
3 x x x x o o o o o o o o o o
4 x x o o o o o o o o o o o o
5 x x o x x x o o o o o o o o
6 x x x x x x x x o o o o o o
7 x o x x x x x x o o o o o o
－－－－－－－－－－－－－－
对ARMA(2,2)过程求eacf表格
x
eacf(x)
运行结果
－－－－－－－－－－－－－－
AR/MA
0 1 2 3 4 5 6 7 8 9 10 11 12 13
0 x x x x x x x x x x x x x x
1 x x o o o x o x o o o o o o
2 x x x o o o o x o o o o o o
3 x x x x o o o o o o o o o o
4 o x x x x o o o o o o o o o
5 x x x x x o o o o o o o o o
6 x x x x o o o o o o o o o o
7 x x x x o x o x o o o o o o
－－－－－－－－－－－－－－
结果分析： 从上面的表格来看，用自相关系数(eacf)来确定ARMA过程的阶数与模型基本相吻合。
4.ARMA模型的参数估计
在确定好ARMA过程的阶数后，可以利用R程序对模型参数进行估计。
下面我们用上证指数的数据来建立ARMA模型。
(1)收集上证指数近99天的数据如下
3190.32 3138.29 3131.11 3136.63 3163.18 3168.90 3160.53 3122.29 3166.65
3133.72 3152.76 3263.48 3280.95 3290.64 3279.25 3269.88 3291.11 3291.38
3310.24 3326.70 3307.17 3288.41 3271.67 3289.64 3256.93 3254.53 3273.75
3259.41 3292.07 3329.57 3289.02 3268.56 3199.16 3184.96 3154.13 3129.85
3262.05 3309.26 3370.65 3487.50 3462.08 3446.98 3480.83 3488.01 3523.00
3558.13 3548.31 3559.47 3546.50 3501.36 3487.86 3474.75 3444.67 3436.59
3410.49 3428.94 3425.34 3421.83 3413.90 3409.48 3391.75 3385.71 3369.11
3348.33 3307.17 3296.38 3275.78 3306.12 3280.46 3297.06 3300.06 3287.61
3296.54 3267.92 3266.14 3292.44 3303.04 3280.81 3322.20 3289.99 3272.05
3293.96 3303.68 3309.62 3317.62 3317.19 3337.86 3333.66 3322.23 3353.82
3351.92 3430.46 3410.50 3392.40 3382.91 3399.25 3402.52 3429.55 3447.84
由于数据是按照时间反序排列，首先做一个变换，把数据按照时间顺序进行排列
x
y
for( i in 1:99) {
y[i]
}
(2)画出时序图然后一阶差分，画出差分后的acf图和pacf图
y
plot(y)
yd
plot(yd)
acf(yd)
pacf(yd)上证指数近99天时序图
一阶差分后的时序图
一阶差分后的acf图
一阶差分后的pacf图
结果分析：显然上证指数的不具有平稳过程的特征。我们进行一阶差分，给出一阶差分后的时间序列acf图和pacf图。
(3)一阶差分后的eacf表
library(TSA)
eacf(yd)
运行结果
-----------------------
AR/MA
0 1 2 3 4 5 6 7 8 9 10 11 12 13
0 o o o x o o o o o o o o o o
1 x o o x o o o o o o o o o o
2 x o o x o o o o o o o o o o
3 x x x x o o o o o o o o o o
4 o x x x o o o o o o o o o o
5 o o o x x o o o o o o o o o
6 x o o x x o o o o o o o o o
7 x o x x o o o o o o o o o o
------------------------------
结果分析：根据eacf表，选取ARMA模型的阶数p=4，q=4.
(4)参数估计与检验建立ARMA(4，4)模型，进行参数估计
out
summary(out)
运行结果：
---------------------------
Call:
arma(x = yd, order = c(4, 4))
Model:
ARMA(4,4)
Residuals:
Min 1Q Median 3Q Max
-115.457 -14.647 3.561 17.975 55.868
Coefficient(s):
Estimate Std. Error t value Pr(>|t|)
ar1 -0.1258 0.2597 -0.484 0.6282
ar2 -0.2730 0.1858 -1.469 0.1417
ar3 -0.1213 0.1853 -0.654 0.5129
ar4 -0.5409 0.2209 -2.449 0.0143 *
ma1 0.2444 0.2817 0.867 0.3857
ma2 0.4421 0.2092 2.114 0.0345 *
ma3 0.4314 0.1956 2.205 0.0274 *
ma4 0.3780 0.2474 1.528 0.1265
intercept -3.9764 7.4762 -0.532 0.5948
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Fit:
sigma^2 estimated as 889, Conditional Sum-of-Squares = 82679.07, AIC = 961.54
－－－－－－－－－－－－－－－－
结果分析：
,
并不显著，其中
的p值最大，进行剔除，然后进行模型估计。这样经过多次系数的筛选，得到以下结果：
R程序实现
out
summary(out)
out
summary(out)
out
summary(out)
out
summary(out)
运行最终结果
--------------------------------------
Call:
arma(x = yd, lag = list(ar = 4, ma = 3))
Model:
ARMA(4,3)
Residuals:
Min 1Q Median 3Q Max
-111.276 -11.932 4.758 18.879 62.097
Coefficient(s):
Estimate Std. Error t value Pr(>|t|)
ar4 -0.31669 0.09945 -3.184 0.00145 **
ma3 0.23203 0.10000 2.320 0.02033 *
intercept -2.60783 3.84574 -0.678 0.49770
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Fit:
sigma^2 estimated as 970.4, Conditional Sum-of-Squares = 90249.36, AIC = 958.13
---------------------------------
结果分析：可以看出，最终模型的系数都显著，且AIC＝958.13，小于最初模型的AIC＝961.54。模型的截距的p值＝0.49770 ，因此可以认为模型截距为0。这样ARMA模型估计为
(5)残差的检验
如果模型拟合充分，则残差序列应该接近白噪声，即均值为0，方差为常数，不存在自相关和偏相关关系。残差检验的工具包括残差时序图、acf、pacf和Ljung－Box检验。
res
plot(res)
acf(res)
pacf(res)残差图
残差的acf图和pacf图
结果分析：从acf图和pacf图可以出，残差是不相关的。但是从时序图中发现个别地方偏差比较大，可以考虑GARCH等模型给予处理。
－－－－－－－－－－－
>Box.test(res,12,"Ljung")#Ljung-Box检验
Box-Ljung test
data: res
X-squared = 10.319, df = 12, p-value = 0.588
－－－－－－－－－－
结果分析：从上面结果可知，p-value = 0.588>0.05，说明残差符合白噪声不相关的假设。
5.ARIMA模型的估计
由于arma函数不能预测，并且考虑到对数据建模时，先对上证指数进行一阶差分，然后进行ARMA模型的参数估计。下面采用arima函数对上证指数建模，结合上面分析的结果，我们确定模型为ARIMA(4,1,3).
(1)模型的估计与检验
R 实现
out1
out1
运行结果为
------------------------------
Call:
arima(x = y, order = c(4, 1, 3))
Coefficients:
ar1 ar2 ar3 ar4 ma1 ma2 ma3
0.1875 -0.1824 -0.0805 -0.2666 -0.0539 0.3221 0.3274
s.e. 0.2551 0.2070 0.2174 0.1184 0.2584 0.1909 0.2115
sigma^2 estimated as 876.9: log likelihood = -471.49, aic = 956.99
----------------------------------
可以看出该模型的aic = 956.99，小于上面的模型AIC＝958.13。
残差的时序图、acf图、pacf图以及Ljung-Box检验
res1
par(mfrow=c(1,3))
plot(res1)
acf(res1)
pacf(res1)
> Box.test(res1,12,"Ljung")
Box-Ljung test
data: res1
X-squared = 4.2965, df = 12, p-value = 0.9775
上面的acf图、pacf图以及Ljung-Box检验说明残差为白噪声，说明模型选取合理。
根据上面的R程序运行的结果可得模型的估计为
(2)预测
下面给出7天的上证指数的预测
pp
print(t(cbind(pp$pred,pp$se)))
运行结果：
-----------------------------
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
pp$pred 3209.124 3224.02774 3234.51432 3218.37592 3207.22367 3203.25815 3203.0522
pp$se 29.612 44.76206 59.00538 74.78575 84.00089 90.61841 96.3358
6.参考文献
【1】R语言与现代统计方法 刘强等编著 清华大学出版社