SAS学习系列39时间序列分析ⅢARIMA模型.docx

1、SAS学习系列39时间序列分析ARIMA模型39. 时间序列分析ARIMA 模型随着对时间序列分析方法的深入研究， 人们发现非平稳序列的确 定性因素分解方法 (如季节模型、趋势模型、移动平均、指数平滑等) 只能提取显著的确定性信息， 对随机性信息浪费严重， 同时也无法对 确定性因素之间的关系进行分析。而非平稳序列随机分析的发展就是为了弥补确定性因素分解方 法的不足。时间序列数据分析的第一步都是要通过有效手段提取序列 中所蕴藏的确定性信息。 Box 和 Jenkins 使用大量的案例分析证明差 分方法是一种非常简便有效的确定性信息的提取方法。而 Gramer 分 解定理则在理论上保证了适当阶数的

2、差分一定可以充分提取确定性 信息。(一) ARMA 模型即自回归移动平均移动模型， 是最常用的拟合平稳时间序列的模 型，分为三类： AR 模型、 MA 模型和 ARMA 模型。一、 AR( p)模型 p阶自回归模型1.模型：xt 0 1xt 1 L pxt p t其中， p 0 ，随机干扰序列 t 为 0 均值、 2 方差的白噪声序列E( t s) 0, ts)，且当期的干扰与过去的序列值无关， 即 E(xtt)=0.由于是平稳序列，可推得均值 0 . 若 0 0 ，称为1 1 L p中心化的 AR(p) 模型，对于非中心化的平稳时间序列，可以令0 (1 1 L p)， xt* xt 转化为中

3、心化。记B为延迟算子， p(B) I 1B L pBp称为 p阶自回归多项式，则 AR(p)模型可表示为： p(B)xt t .2.格林函数用来描述系统记忆扰动程度的函数， 反映了影响效应衰减的快慢 程度(回到平衡位置的速度) ，Gj 表示扰动 t-j 对系统现在行为影响的 权数。例如， AR(1)模型(一阶非齐次差分方程) ，Gj 1j, j 0,1,2,L模型解为 xt Gj t j .j03.模型的方差2 对于 AR(1)模型， Var (xt ) G j2Var ( t j) 2 .j 0 1 14.模型的自协方差 对中心化的平稳模型，可推得自协方差函数的递推公式：用格林函数显示表示：

4、i0j0对于 AR(1) 模型，5.模型的自相关函数递推公式：对于 AR(1)模型， (k) 1k (0) 1k .平稳 AR(p)模型的自相关函数有两个显著的性质：1)拖尾性指自相关函数 (k始) 终有非零取值， 不会在 k 大于某个常数之后 就恒等于零；2)负指数衰减随着时间的推移， 自相关函数 (k会) 迅速衰减，且以负指数 (ik 其 中 i 为自相关函数差分方程的特征根)的速度在减小。6.模型的偏自相关函数自相关函数 （k实） 际上并不只是 xt与 xt-k 之间的相关关系， 它还会 受到中间 k-1个随机变量 xt-1, , xt-k+1的影响。为了能剔除了中间 k-1 个随机变量

5、的干扰，单纯测度 xt与 xt-k之间的相关关系，引入了滞后 k偏自相关函数（ PACF），计算公式为：其中，滞后 k偏自相关函数实际上等于 k阶自回归模型第 k个回归系数 kk ：两边同乘以 xt-k，求期望再除以 （0） 得到取前 k 个方程构成的方程组：称为 Yule-Walker 方程，可以解出 kk可以证明平稳 AR(p)模型，当 kp 时， kk 0. 即平稳 AR(p)模 型的偏自相关函数具有 p 步截尾性。注：实际上样本的随机性使得偏自相关函数不是严格截尾， 例如 上面两图都 1 阶显著不为 0， 1 阶之后都近似为 0.二、MA(q) 模型 q阶移动平均模型1.模型：其中，

6、q 0 ，随机干扰序列 t 为 0 均值、 2方差的白噪声序列 (E( t s) 0, ts)。若 =0，称为中心化的 MA(q) 模型，非中心化的 MA(q) 模型可以 通过 xt* xt 转化为中心化。记 B 为延迟算子， q(B) I 1B L qBq称为 q 阶自移动平 均系数多项式，则中心化 MA(q)模型可以表示为 xt q(B) t .2.模型的方差3.模型的自协方差 只与滞后阶数 k 相关，且 q阶截尾。当 k=0 时，当 1k q 时，当 kq 时， (k) 0.4.模型的自相关函数： (k) (k)(q 阶截尾性)(0)5.模型的滞后 k 阶偏自相关函数(中心化)可以证明滞

7、后 k 阶偏自相关函数具有拖尾性。6.模型的可逆性以 MR(1) 为例，模型： xt t 1 t 1 或 x1t t1 1 1 B1它们的自相关函数 1 1 /(1 12 )相同(即相同的自相关函数对应不 同的回归模型)，为了保证对应的唯一性，需要增加约束条件，即 MR(q) 模型的可逆性条件。观察两个模型的第二种表示：当 | 1 | 1时，模型收敛、模型 不收敛；当 | 1 | 1时，模型不收敛、模型收敛。表示成收敛形式的 MR(q)模型称为可逆 MR(q) 模型。一个自相关 函数只对应唯一一个可逆 MR(q) 模型。三、 ARMA(p, q) 模型自回归移动平均模型1.模型其中， p 0，

8、 q 0，随机干扰序列 t为 0 均值、 2方差的白噪声 序列( E( t s) 0 , ts)，且当期的干扰与过去的序列值无关，即 E(xtt)=0.若 0 =0 ，则称为中心化的 ARMA(p,q) 模型。引入延迟算子，中 心化的 ARMA(p,q) 模型可表示为： p(B)xt q(B) t .显然， AR(p)和MA(q)模型是 ARMA(p,q)模型的特例。2.数字特征1)均值： E(xt )2)自协方差函数： (k) 2 GiGi k ，其中 Gi 为格林函数；i0(k) GiGi k3)自相关函数： (k) (k) i 0(0) 2i03. 模型的初步定阶对于平稳非白噪声序列，计

9、算出样本自相关系数( ACF )和偏自 相关系数(PACF)，根据其性质估计自相关阶数 p? 和移动平均阶数 q?， 称为 ARMA(p,q) 模型的定阶。可以推导出：样本自相关函数 ?(k) 和偏自相关函数 ?kk都近似服1从正态分布 N(0, ).n取显著水平 =0.05，若样本自相关系数和样本偏自相关系数在最 初的 k 阶明显大于 2 倍标准差，而后几乎 95%的系数都落在 2 倍标准 差的范围内， 且非零系数衰减为小值波动的过程非常突然， 通常视为 k 阶截尾；若有超过 5%的样本相关系数大于 2 倍标准差，或者非零 系数衰减为小值波动的过程比较缓慢或连续，通常视为拖尾。4.参数估计对

10、非中心化的 ARMA(p,q) 模型xtq(B) p(B)参数 可用样本均值来估计总体均值（矩估计法） ，初步定阶估计出 自相关阶数 p? 和移动平均阶数 q?后，模型共有 p+q+1 个未知参数： ,L , , ,L , , 2 .1 p 1 q（1）参数的矩估计用时间序列样本数据计算出延迟 1阶到 p+q 阶的样本自相关函数 ?（k） ，延迟k阶的总体自相关函数为 k（ 1,L , p, 1,L , q）. 用计算出 的样本自相关函数来估计总体自相函数，得到 p+q 个联立方程组： 从中解出 1,L , p, 1,L , q的值作为未知参数估计值 ?1,L , ?p , ?1,L , ?q

11、. ARMA（p,q） 模型的两边同时求方差，并把前面的参数的估计值代入， 可得白噪声序列的方差估计为：（2）参数的极大似然估计当总体分布类型已知时， 极大似然估计是常用的估计方法。 其基 本思想是，认为样本来自使该样本出现概率最大的总体。因此，未知参数的极大似然估计，就是使得似然函数（即联合密 度函数）达到最大值的参数值：在时间序列分析中， 序列的总体分布通常是未知的。 为了便于分 析和计算， 通常假设序列服从多元正态分布， 它的联合密度函数是可 导的。在求极大似然估计时，为了求导方便，常对似然函数取对数， 然后对对数似然函数中的未知参数求偏导数， 得到似然方程组。 理论 上，只要求解似然方

12、程组即可得到未知参数的极大似然估计。 但在实 际上是使用计算机经过复杂的迭代算法求出未知参数的极大似然估 计。两种估计的比较：矩估计的优点是不要求知道总体的分布， 计算量小， 估计思想简 单直观。但缺点是只用到了样本自相关系数的信息， 序列中的其他信 息被忽略了，这导致估计精度一般较差。因此，它常被作为极大似然 估计和最小二乘估计的迭代计算的初始值。极大似然估计的优点是充分应用了每一个观察值所提供的信息， 因而它的估计精度高，同时，还具有估计的一致性、渐近正态性和渐 近有效性等优良统计性质，是一种非常优良的参数估计方法。( 3)参数的最小二乘估计使 ARMA(p,q) 模型的残差平方和达到最小

13、的那组参数值： 通过计算机借助迭代方法求出。 由于充分利用了序列的信息， 该方法 估计精度最高。在实际运用中， 最常用的是条件最小二乘估计， 假定时间序列过去未观察到序列值等于序列均值，可得到残差的有限项表达式：于是残差平方和达到最小的那组参数值为：5.模型和参数的显著性检验ARMA(p,q) 模型中，使用 QLB 统计量检验残差序列的自相关性， 为了克服 DW 检验的有偏性， Durbin 在 1970 年提出了修正的 Durbin h 统计量：其中，n 为观察值序列的长度， 2 为延迟因变量系数的最小二乘估计 的方差。参数的显著性检验是要检验每一个模型参数是否显著非零。 若某 个参数为零，

14、模型中包含这个参数的乘积项就为零，可以简化模型。 因此，该检验的是为了精简模型。原假设 H0：某未知参数 j=0；H1：j0. 可以构造出检验未知 参数显著性的 t(n-m)检验统计量，其中 m 为参数的个数。6.模型优化当一个拟合模型在置信水平 下通过了检验， 说明了在该置信水 平下该拟合模型能有效地拟合时间序列观察值的波动。 但是这种有效 的拟合模型并不是惟一的。如果同一个时间序列可以构造两个拟合模型， 且两个模型都显著 有效，那么应该选择哪个拟合模型用于统计推断呢？通常采用 AIC 和 SBC 信息准则来进行模型优化。（1）AIC 准则最小信息量准则由日本统计学家赤池弘次（ Akaike

15、）于 1973 年提出，是一种考 评综合最优配置的指标，它是拟合精度和参数未知个数的加权函数：AIC =2ln（ 模型中极大似然函数值 ）+2（ 模型中未知参数个数 ） 使其达到最小值的模型被认为是最优模型。（ 2）BIC/SBC 准则AIC 准则的不足：若时间序列很长，相关信息就越分散，需要多 自变量复杂拟合模型才能使拟合精度比较高。在 AIC 准则中拟合误 差等于nln（ ?2） ，即随样本容量 n 增大，但模型参数个数的惩罚因子 （始终=2）却与 n无关。因此在样本容量 n 趋于无穷大时，由 AIC 准则选择的拟合模型不收敛于真实模型， 它通常比真实模型所含的未 知参数个数要多。为了弥补

16、 AIC 准则的不足， Akaike于1976年提出 BIC准则。而 Schwartz在1978年根据贝叶斯理论也得出同样的判别准则， 称为 SBC 准则。 SBC 准则定义为：SBC=2ln（模型中极大似然函数值 ）+ln（ n）（模型中未知参数个数 ） 即将未知参数个数的惩罚权重由常数 2 变成了 ln（n）。在所有通过检 验的模型中使得 AIC 或 SBC 函数达到最小的模型为相对最优模型 （因 为不可能比较所有模型） 。7.模型预测 即利用时间序列已观察到的样本值对时间序列在未来某个时刻 的取值进行估计。常用的预测方法是线性最小方差预测。根据 ARMA(p,q) 模型的平稳性和可逆性，

17、可以用格林函数的传 递形式和逆转函数的逆转形式等价描述该序列： 右式代入左式得：xt Gi I j xt i j Gi I j xt i j Ci xt 1 ii 0 j 0 i 0 j 0 i 0可见， xt是历史数据 xt-1, xt-2, 的线性函数。对于任意一个将来时刻 t+l，也可以用上式预测，但 xt+l -1, x, t+1 未知。根据线性函数的可加性， 所有未知信息都可以用已知信息的线 性函数表示出来，并用该线性函数进行估计： 用et(l) xt l x?t l 来衡量预测误差，最常用的预测原则是预测误差的 方差最小法：在线性预测方差最小法下得到的估计值 x?t l 是在序列

18、xt, xt-1, 已 知的情况下得到的条件无偏最小方差估计值。 且预测方差只与预测步 长l 有关，而与预测起始点 t无关。预测步长 l 越大预测值的方差越大，因此只适合于短期预测。在 正态假定下，估计值 x?t l的 1-的置信区间为：二) ARIMA 模型混和自回归移动平均模型、原理也称 Box-Jenkins 模型，用来处理单变量同方差的非平稳时间序ARIMA(p,d,q) 模型的形式如下：其中， d (I B)d为 d阶差分，为平稳可逆 ARMA(p,q) 模型的自回归和移动平均系数多项式。可见， ARIMA 模型的实质就是差分运算与 ARMA 模型的组合 任何非平稳序列只要通过适当阶

19、数的差分实现平稳， 就可以对差分后 序列进行 ARMA 模型的拟合了。d 阶差分后的序列可表示为：其中， Cdi 为组合数，即 d 阶差分后序列等于原来序列的若干序列值 的某种加权和。、建模步骤分为三个阶段：识别阶段、估计阶段和预测阶段1. 识别阶段使用 identify 语句来指定响应变量序列并且识别候选 ARIMA 模 型。一般先对序列进行非线性、差分和平稳性检验，可能对序列进行 差分，然后计算自相关系数 ACF、逆自相关系数 IACF 、偏自相关系 数 PACF 和互相关系数。 此阶段的输出通常会建议一个或多个可拟合 的 ARIMA 模型。如果模型确定，还可以检验样本自相关系数 SACF

20、 和样本偏自相关系数 SPACF，以分出模型的类型。2. 估计阶段使用 estimate 语句来指定 ARIMA 模型去拟合在前面 identify 语 句中指定的响应变量，并且估计该模型的参数。 estimate 语句也生成 诊断统计量从而帮助判断该模型的适用性。关于参数估计值的显著性检验可以指出模型里的一些项是否不 需要：拟合优度统计量 R2 可帮助比较该模型和其他模型的优劣；白 噪声残差检验可指明残差序列是否包含可被其他更复杂模型采用的 额外信息， 如果诊断检验表明模型不适用， 则可尝试另一个模型然后 重复估计和诊断。3.预测阶段使用 forecast 语句来预测时间序列的未来值，并对这

21、些来自前面 estimate 语句生成的 ARIMA 模型的预测值产生置信区间。(3)PROC ARIMA 过程ARIMA 过程采用 Box-Jenkins 方法建立模型，是集一元时间序列 模型判定、参数估计和预测为一体的多功能综合工具。当 ARIMA 模 型包括其他时间序列作为输入变量时， 有时也被称为 ARIMAX 模型。 ARIMA 模型还支持干预或中断时间序列模型、 误差的多元回归分析、 任意复杂程度的有理转移函数模型。基本语法：proc arima data=数据集 out=输出数据集 ;where 条件表达式 ;identify var=变量 () ;estimate ;forec

22、ast ;说明：(1)where 语句指定用于分析的时间间隔， 通常条件表达式是有关日期变量的条 件表达式，例如： 31dec98d日 期变量 31dec99d( 2)identify 语句主要完成时间序列的差分计算，样本 ACF、IACF 和 PACF 函数 的计算、卡方检验统计量和白噪声自相关检验的 p 值的计算。 主要选 项：var=变量(d1, ,dk) 是必选项，指定要分析的时间序列变量， 按括号内列出的差分周期列表来计算时间序列的滞后差分。例如：var=X(1) 为对滞后 1 项的序列差分，即 X t-X t- 1;var=X(2) 为对滞后 2 项的序列差分，即 X t-X t-

23、 2; var=X(1,1) 为 X 进行二阶差分，即 (Xt-Xt-1)-(Xt-1-Xt-2); nlag=数字指定计算自相关的滞后数， 其值应大于 p+d+q，小 于观测数，默认值为 24；crosscorr=(干预变量 (d1)列出有 var=指定的响应序列的交叉 相关变量。 干预变量在交叉相关变量中。 交叉相关变量的差分由圆括 号内的差分滞后数确定。( 3)estimate语句对已执行的 identify 语句中的响应变量规定一个模型， 主要选项： p=(p1,p2, )(p1,p2, )定义一个在 p 中指定的滞后处具有自回归参数的模型， p 的默认值为 0;q=(q1,q2, )

24、(q1,q2, )定义一个在 q 中指定的滞后处具有 滑动平均参数的模型， q的默认值为 0。如果 p=和q=都没有指定，则 拟合随机模型；noconstant 在模型中舍弃常数项 ;noint 在该模型中不拟合截距参数；method=ml | uls | cls 指定估计时使用的方法，分别为极大似然方法、无条件最小二乘法、有条件最小二乘法，默认为 cls;outest=数据集将参数估计值输出到指定的数据集； outmode=数据集将模型和参数估计值输出到指定的数据集； outstat=数据集将模型诊断统计量输出到指定的数据集； plot 可以绘制残差自回归函数等；(4)forecast语句利

25、用 estimate语句所产生的参数估计生成时间序列的预测值， 主要选项：alpha= 设置预测置信限的大小，上下置信限的置信水平为1-，默认值为 0.05;lead=n 指定要计算的多步向前预测值的次数，默认为 24； back=n指定在数据末尾前 n 个观测值开始进行分步预测，默 认为 0;interval=时间间隔指定观测之间的时间间隔， 常用的时间间隔为 year、qtr、 month、week、weekday、day、hour、minute、second; id= 变量名指明输入数据集中一个变量， 用于识别与观测有关的时间周期；out=数据集将预测值和其他值输出到一个指定的数据集中例

26、 1 有一组 1949 年至 1961 年国际航线旅客月度人数的记录：YEAR123456789101112194911211813212912113514814813611910411819501151261411351251491701701581331141401951145150178163172178199199184162146166195217118019318118321823024220919117219419531961962362352292432642722372111802011954204188235227234264302293259229203229195524

27、223326726927031536434731227423727819562842773173133183744134053553062713061957315301356348355422465467404347305336195834031836234836343549150540435931033719593603424063964204725485594634073624051960417391419461472535622606408461390432使用 ARIMA 过程进行建模和预测一) 读入数据、绘制时间序列图、判断平稳性， 代码：data arimad01;date=int

28、nx( month , 31dec1948d ,_n_);input x ;format date monyy5. ; datalines ;112118132129121135148148136119104118115126141135125149170170158133114140145150178163172178199199184162146166171180193181183218230242209191172194196196236235229243264272237211180201204188235227234264302293259229 203 229242233267269270315364347312274237278284277317313318374413405355306271306315301356348355422465467404347305336340318362348363435491505404359310337360342406396420472548559463407362405417391419461472535622606408461 390 432run ;procsgplotdata =arimad01;series