时间序列分析-第六章 ARMA模型的参数估计

1、第六章第六章 ARMAARMA模型的参数估计模型的参数估计n第一节 AR(p)模型的参数估计n第二节 MA(q)模型的参数估计n第三节 ARMA(p,q)模型的参数估计n第四节 求和模型及季节模型的参数估计 第一节第一节. AR(p). AR(p)模型的参数估计模型的参数估计n目的：为观测数据建立AR(p)模型 (1.1) 假定自回归阶数p已知，考虑回归系数 和零均值白噪声 的方差 的估计。n数据 的预处理：如果样本均值不为零，需将它们中心化，即将它们都同时减去其样本均值 再对序列按(1.1)式的拟合方法进行拟合。 tptptttXXXX2211Tp),(1t2nxxx,21nttnxnx1/

2、1 n假定数据 适合于以下模型 (1.2)其中，p为给定的非负整数， 为未知参数，记 为系数参数， 为独立同分布序列，且 ， 与 独立，参数满足平稳性条件。nxxx,21nptXXXXtptpttt, 1,2211p,21Tp),(1t422, 0tttEEEt,tsxs A. AR(p)A. AR(p)模型参数的模型参数的Yule-WalkerYule-Walker估计估计n对于AR(p)模型，自回归系数 由AR(p)序列的自协方差函数 通过Yule-Walker方程 唯一决定，白噪声方差 由 决定。ppppppaaarrrrrrrrrrrr2102120111021prrr,102jpjj

3、rr102 nAR(p)模型的自回归系数和白噪声方差的矩估计 就由样本Yule-Walker方程 (1.3) 和 (1.4) 决定。21,),(Tppppppprrrrrrrrrrrr2102120111021jpjjrr102 n令 则(1.3),(1.4)式可写为ppppppppprrrrrrrrrrrr,2121021201110bpppb n实际应用中，对于较大的p,为了加快计算速度可采用如下的Levison递推方法 递推最后得到矩估计pkkjaaaaarrarrararjkkkkjkjkkjkjkjjkjjkkkkkkkk,1)()1 (/1,1, 1, 11110111, 12,2

4、122011102022,2,1 ,1,),(),(pTppppTpaaa上式是由求偏相关函数的公式：导出。kkkkkkkkkaaa2121211121111 n定理1.1 如果AR(p)模型中的 是独立同分布的 ，则当 时(1) (2) 依分布收敛到p维正态分 布 。 t42), 0(tEWNn22,pjpjTppn),(11),(12pN0 n注：用 表示 的第 元素时，可知 依分布收敛到 ，于是 的95%的渐近置信区间是 在实际问题中， 未知，可用 的 元素 代替 ，得到 的近似置信区间jj,12pjj)(jjn), 0(, jjNj/96. 1,/96. 1,nnjjjjjjjj,12

5、pjjjj,jj,j/96. 1,/96. 1,nnjjjjjj B. AR(p)B. AR(p)模型参数的最小二乘估计模型参数的最小二乘估计n如果 是自回归系数 的估计，白噪声 的估计定义为 通常 为残差。n我们把能使 (1.6) 达到极小值的 称为 的最小二乘估计。p,21p,21jnjpxxxxpjpjjjj12211),(njpj1,npjptptttxxxxs122211)( n记 则 ，于是 的最小二乘估计为 即 ,21211121pnnnppppnppxxxxxxxxxxxxxyyyxyyxxxTTTTTTs)(yxxxTT1)()(inf)()(1yxxxxyyyssTTTT

6、n相应地，白噪声方差 的最小二乘估计 式中 为 的p个分量。nptptpttTTTTxxxpnpnspn121112)(1)(1)(1yxxxxyyy2p,21 n定理1.2 设AR(p)模型中的白噪声 是独立同分布的， 是自回归系数 的最小二乘估计，则当 时，依分布收敛到p维正态分布 n注：对于较大的n，最小二乘估计和矩估计(Yule-Walker)估计的差别不大。t),( ,214ptEp,21n),(2211ppn),(12pN0 .),/1(模型，有则对为最小二乘估计，估计，ker为记).(就称),1(/如果).1(是依概率有界的，记为就称,)|(|sup，使得存在正数,0任何是非零常

7、数列，如果对是时间序列，：设1.1定义nnOARWalYulecOOcOMPMcpLLnpnpnnpnnnnnn C. AR(P)C. AR(P)模型的极大似然估计模型的极大似然估计n假定模型AR(p)中的 为正态分布，则观测向量 的高斯似然函数为 相应的对数似然函数为 其中， 为 的协方差阵， 表示 的行列式，使得对数似然函数达到极大值的 和 称为 和 的极大似然估计。tTnnxxx),(x21 )xx21exp(|)2(),|,(1212212nnTnnnnxxxL nnTnnnnxxxlxx21|21)2log(2),|,(121212 |nnTnxxx),(21n),|,(212nxx

8、xl 22n从另一角度考虑：.)2ln(2,)(21)ln(2)(21)ln(2),(.)(21exp)2(),(,).21exp()2(,22121222121222112221是是常常数数其其中中为为相相应应的的对对数数似似然然可可定定义义的的似似然然函函数数于于是是可可得得基基于于有有联联合合密密度度函函数数服服从从正正态态分分布布，则则由由于于 pNccSpNcxxpNlxxLxxnptpjjtjtnptpjjtjtpnnnpttpnnpt 的的最最小小二二乘乘估估计计。的的最最小小值值点点，从从而而是是的的最最大大值值点点实实际际上上是是容容易易看看出出，是是常常数数这这里里，表表达

9、达式式，得得到到将将上上式式代代入入于于是是，得得的的最最大大值值点点，解解方方程程为为求求)(),(.)(21)(ln2),(),().(10)(212),(),(200222242222SlccSSpNllSpnSpnll n注：当n充分大时，AR(p)模型参数的极大似然估计、最小二乘估计和矩估计(Yule-Walker估计)三者都非常接近，即三者渐近相等，它们都可以作为AR(p)模型的参数估计，这是AR(p)模型的独有的优点。 n例1.1. 由下列AR(1)序列产生长度为n=300的样本，计算出前5个样本自协方差函数值为求参数的矩估计和最小二乘估计。(1) 参数 的矩估计 分别为将样本自

10、协方差函数值代入得) 1 , 0(,5 . 01NXXtttt0123. 0,1773. 0,3886. 0,7771. 0,5419. 143210rrrrr21,21,11021011,/rrrr150. 1,504. 021 (2) 参数 的最小二乘估计 分别为21,21,074. 1,506. 021 n例1.2 求AR(2)模型 参数 的估计，这里n=300, (1) AR(2)模型的矩估计为ttttXXX2211221,22110221202120221202011)(rrrrrrrrrrrrr,24. 0, 121) 1, 0(. .Ndi it 计算出的前5个样本协方差函数值为

11、将其值代入上式得：(2) 最小二乘估计2705. 0,8060. 0,4362. 1,2171. 2,7888. 243210rrrrr922433. 0,318064. 0,047842. 1221939530. 0,328336. 0,071838. 1221 n注：一般在求高阶AR(p)模型参数的矩估计时，为了避免求高阶逆矩阵，可采用求偏相关函数的递推算法，求出 即为 的矩估计，将它们代入 的表达式可得 。ppppaaa,21p,2122 D. AR(p)D. AR(p)模型的定阶模型的定阶1. 1. 偏相关函数的分析方法偏相关函数的分析方法n一个平稳序列是AR(p)序列当且仅当它的偏相

12、关函数是p步截尾的。n如果 p步截尾：当 时， ； 而 ，就以 作为p的估计。,kkapk0,kka0, ppap n定理1.3 设 由 定义，如果AR(p)模型中的白噪声是独立同分布的， ，则对确定的kp,当 时， 依分布收敛到k维正态分布 。kjaaaaaaarrarrarrajkkkkjkjkkjkjkjjkjjkkkjkjkjjkjjkkkk,2 , 1)1)()(1,1, 1, 1111111110111, 110111kkkkaaa,2,1 ,4tEn),(,2,2,1 ,1 ,kkkkkkkkaaaaaan),(12kN0 n推论：在定理1.3的条件下，对kp, 依分布收敛到标准

13、正态分布N(0,1)。n根据推论，对于AR(p)序列和kp，当样本量n比较大时， 以近似于0.95的概率落在区间 之内。于是对于某个固定的k，以 作为p的估计。 kkan,kkan,/96. 1 ,/96. 1nn1,/96. 1|:|sup,kjnajpjj n或者根据推论有如下的检验方法：对于某个正整数p, 显著地异于零，而 近似等于零,其满足 (或 )的个数占 的比例近似地为68.3%(或95.5%),则近似地认为 在p步截尾， 初步判定为AR(p)。 ppa,00,2,21, 1,ppppppaaanakk/1|,nakk/2|,pp 0,kkatX n例1.3（例1.1续）使用样本偏

14、相关函数对AR(p)的模型阶数作初步的判定。结果：取上限 ，样本自相关函数 呈拖尾状，而从15个偏相关函数来看,除 显著异于零之外，其余14个中绝对值不大于 的有10个，于是结论：初步判定为AR(1)模型。150pk0577. 0300/1/1n%3 .68%43.7114/101 , 1 a n前15个样本偏相关函数-0.20.02468101214 n例1.4（例1.2续）使用样本偏相关函数对AR(p)的模型阶数作初步的判定。结果：取上限 ，样本自相关函数 呈拖尾状，而从15个偏相关函数来看，除 显著异于零之外，其余14个中绝对值不大于 的有9个,于是结论：初步判定为AR

15、(2)模型。2,21 , 1,aa150pk0577. 0300/1/1n%3 .68%2 .6913/9 n前15个样本偏相关函数-0.4-0.20.00.81.02468101214 2. 2. AICAIC准则方法准则方法( (A-Information Criterion)A-Information Criterion)n为了使拟合残差平方和 尽量小，而又不至于引入过多的虚假参数的估计，Akaike于1973年引入如下的准则函数，假定已有阶数p的上阶 ， AIC(k)的最小值点 (若不唯一，应取小的)称为AR(p)模型的AIC定阶，即)(s0P02, 1 , 0,2)

16、(log)(PknkkkAICp )(min)(00kAICpAICPk n具体步骤：1. 取定p=k时，根据数据 使用前一小节所提的任何一种参数的估计方法，给出噪声方差 的估计 ；2. 再找出AIC取极小值时，所对应的阶数p.n注：AIC定阶并不相合，AIC定阶通常会对阶数略有高估。故在应用中，当样本量不是很大时，使用AIC定阶方法。nxxx,2122 n为了克服AIC定阶的不相合性，可使用BIC准则方法。设 为AR序列，则BIC准则函数为 将此准则函数达到最小值的解 作为p的估计，就是BIC准则方法。n注： 1. 理论上已证明BIC准则的定阶具有相合性。 2. 当n不是很大时，用BIC定阶

17、有时会低估阶数p，造成模型的较大失真，故在实际问题中，特别当样本量不是很大时,BIC的定阶效果并不如AIC定阶准则。tX02, 1 , 0,log)(log)(PknnkkkAICp n例1.5（例1.1续）n=300个观测，定阶。方法：观察偏相关函数，确定上界是P=10，对p=1,2,10分别解Yule-Walker方程得到 的Yuler-Walker估计，再对p=1,2,10分别计算出AIC和BIC函数，计算结果如下： p12345AIC(p)2.98392.99393.00383.00023.0050BIC(p)2.99963.02523.05083.06293.08342 结果：AIC

18、(1)和BIC(1)分别是AIC和BIC函数的最小值。 结论：由AIC和BIC定阶可知阶数p=1.p678910AIC(P)3.01103.01733.02663.02833.0308BIC(P)3.10513.12713.15203.16943.1876 nAIC函数图 2.982.993.003.013.023.033.0412345678910AICX1 nBIC函数图2.953.003.0012345678910BICX1 n例1.6 （例1.2续） n=300个观测，定阶。方法：观察偏相关函数，确定上界是P=10,对p=1,2,10分别求出 的估计，再对p=1

19、,2,10，计算AIC和BIC函数，计算结果如下：p12345AIC(p)2.84702.72772.73682.72812.7377BIC(p)2.86272.75912.78382.79382.81612 结果：AIC(2)和BIC(2)分别是AIC和BIC函数的最小值。结论：由AIC和BIC定阶可知阶数p=2。 p678910AIC(p)2.74652.75672.75922.76272.7688BIC(p)2.84062.86652.88462.90382.9256 nAICAIC函数图函数图 2.722.742.762.782.802.822.842.8612345678910AIC

20、X2 nBIC函数图2.752.802.852.902.9512345678910BICX2 n例1.7：独立重复1000次实验，每次产生符合模型AR(4) 的300个观测，得到AIC和BIC定阶情况如下：12345678910AIC定阶052256741136129211411BIC定阶145559476720000, 2 , 1,18. 011. 037. 016. 14321tXXXXXtttttt n在1000次模拟计算中AIC将阶数定为4的有674次，而BIC阶数定为4的有476次。BIC定阶对阶数低估的比率为51.5%n增大样本量n=1000，获得如下结果：12345678910A

21、IC定阶0007391244537251218BIC定阶041990500000 nAIC定出的平均阶数是Avc(AIC)=4.593，BIC定出的平均阶数是Avc(BIC)=3.996，故对于较大的样本量有必要综合考虑AIC定阶和BIC定阶。 E. E. 拟合模型的检验拟合模型的检验 现有数据 ，欲判断它们是否符合以下模型式中 被假定为独立序列,且 与 独立。n原假设 ：数据 符合AR(p)。故在 成立时，下列序列 为独立序列 的一段样本值序列。nxxx,212, 1,2211pptXXXXtptpttt422, 0tttEEEtt,tsXs0Hnxxx,210HnptXXXXptptttt

22、, 1,2211t n步骤：1. 首先，根据公式 计算出残差的样本自相关函数，2. 利用上一章关于独立序列的判别方法，判断 是否为独立序列的样本值3. 根据判断结果，如果接受它们为独立序列的样本值，则接受原假设，即接受 符合AR(p),否则，应当考虑采用新的模型拟合原始数据序列。2 , 1 , 0),(/)()(2 , 1 , 0,1)(01krrkpnrkkkpntkptptknp,1nxxx,21 n例1.8（例1.5续） 拟合后，给出残差头15个数据，有11个落在之间， 故不能否定原假设，即 符合AR(1)模型。299/1 ,299/1%3 .68%33.7315/11nxxx,21 n

23、残差的图形-4-202450100150200250300RESIDX1 n残差的自相关函数-0.15-0.10-0.050.000.050.102468101214 n例1.9（例1.6续） 拟合后，给出残差头15个数据，有15个落在 之间，故不能否定原假设，即 符合AR(2)模型。298/96. 1 ,298/96. 1nxxx,21 n残差的图形-3-2-1012350100150200250300RESIDX2 n残差的自相关函数-0.15-0.10-0.050.000.050.100.152468101214 第二节第二节 滑动平均模型拟合滑动平均模型拟合n对于已给的时间序列数据 ，

24、用MA(q)式的滑动平均模型去拟合它们，称为滑动平均模型拟合。n滑动平均模型拟合主要包括：(1) 判断滑动平均模型MA的阶数；(2) 估计模型的参数；(3) 对拟合模型进行检验。nxxx,21 一一. . 参数估计参数估计n假定数据序列 适合以下模型 (2.1)其中 为独立同分布的序列,且 ,q为给定的非负整数， 为未知参数，并满足可逆性条件。nxxx,21ntXqtqtttt,2 , 1,2211t22, 0ttEE4tETq),(21 1. 1. 参数的矩估计方法参数的矩估计方法nMA(q)序列的自协方差函数与MA(q)的模型参数有如下公式： 故， 和 的矩估计 和 ，为 (2.2)kqi

25、kkiikqjjqkrr1201221)()1 (当22kqikkiikqjjqkrr1201221)()1 (当 (1) (1) 解析法解析法n对于阶数较低的MA(q)模型，例如MA(1)和MA(2)，可利用解析法求解。n对于MA(1)模型： , 和 满足 可得 和 的矩估计分别为11tttX122112210)1 (rr12)411 (2,411221022111r n例4.11 由MA(1)模型产生长度n=300的样本，计算出前两个样本自协方差函数值 ,由上述讨论) 1, 0(,5 . 01NXtttt487. 0,220. 110rr9776. 0,4979. 021 n对于MA(2)

26、模型： , 其中 满足 可得 的估计为:当 时 221,2211ttttX22221121222120)()1 (rrr221,021)21(214121)21(21412122122222122222 当 时，从而可得，021)21(214121)21(21412122122222122222220222121,)1 (r n例4.12 求MA(2)模型的n=950的样本的参数 的矩估计。解：已知前三项的样本自相关函数分别为使用上述公式，可得到如下估计值) 1 , 0(,3 . 05 . 021NXttttt221,246. 0,221. 0,308. 1210r040. 1,4025. 0

27、,3094. 0212 (2). (2). 线性迭代算法线性迭代算法n将(2.2) 式表示为 (2.3)在可逆域内，给定 的初值,代入(2.3)式右边，得到一步迭代值 ，再将它们代入(2.3)式右边，得出(2.3)式左边的第二不迭代值 ，同法重复直到某步 ，kqikiiqjjkkqjjqkrrr112012021 ,)1 (1221),(和Tq) 1 () 1 (2和)2()2(2和) 1() 1(2mm和 设有精度 ,当同时成立时，就停止迭代（否则继续迭代下去），以 作为 的矩估计。21,2221| )() 1(|, 2 , 1,| )() 1(|mmqjmmjj)()(2mm和2和 (3)

28、 (3) Newton- Newton- RaphsonRaphson 算法算法n优点：方法简便、收敛速度快n缺点：使用该算法得到的解 不能保证满足属于可逆域，需要采用调整方法才可做到。详见时间序列的分析与应用或应用时间序列分析。q,21 2. 2. 极大似然估计极大似然估计n若(2.1)中， 为正态分布，则 服从 分布，其中 是 的协方差矩阵。于是有似然函数：其中， 。使似然函数达到极大值之解的 和 ，即为 和 的极大似然估计。tnxxx,21),(nN0nnxxx,21nnTnnnnxxxlxx12121221|21)2log(2),|,(Tnnxxx),(21x22 n近似极大似然估计方

29、法：假定(2.1)式中的初值 给定，不妨设为零值。则由(2.1)式和数据 可以求出 (2.4) 于是可得到如下近似似然函数为： (2.5) 记 q110,nxxx,21nkxqkqkkkk,2 , 1,2211nkknnl122222122loglog),(nkkS12)( 由(2.5)式决定的近似极大似然估计 和 满足以下方程于是 为以下方程的解而,)(,)(,00222222422SlSnnl2, 0)(SnS)(2 3. 3. 自回归逼近方法自回归逼近方法n原理：可逆的MA(q)模型有逆转形式 模型，且逆转形式中的无穷阶自回归系数满足以指数衰减到零的趋势，故一个可逆的MA模型可用适当高阶

30、的AR模型近似。n用一个高阶的AR模型拟合一个较低的MA序列称为自回归逼近拟合方法。)(AR n步骤：1. 对原始数据 进行自回归模型拟合。可用AIC定阶，求参数的Yule-Walker估计，在进行检验；或直接拟合AR(p)模型。其中，当n不太大时，取 ；当n很大时，取 。将拟合后模型记为 (2.6)nxxx,2110/nP ncPlogtptptttxxxx2211 2. 利用(2.6)式，计算拟合残差：于是(2.1)式的模型可近似写为 记 (2.7)nPtxxxxptptttt, 1,2211nqPtXqtqtttt, 1,2211nqPqPqnnnPqPqPPqPqPnqPqPExxx2

31、121211121,x 于是，(2.7)可简记为故， 和 的最小二乘估计分别为和n优点：不涉及非线性代数方程，易于实际应用。xE2xTTEEE)(1nS/)(2 二二. . 阶数的估计阶数的估计1. 1. 自相关函数估计方法自相关函数估计方法n依据：一个平稳序列为MA(q)序列的充要条件是它的自协方差函数q步截尾。对于MA(q)模型，当kq，n充分大时， 的分布渐近正态 ，于是当kq，n充分大时，下列等式近似成立)21(1,0(12 qjjnN %5.95)212|(|%3.68)211|(|1212qjjkqjjknPnP或k n方法：对于每一个正整数q,计算样本自相关函数(M一般取为 左右

32、),考察其中满足 的个数是否占M的68.3%(或95.5%)左右,如取某 显著地异于零，而 近似等于零，并满足上述不等式的个数达到了68.3%(或95.5%)左右比例，则初步认为 在 步截尾， 初步判定为n)212|(211|1212qjjkqjjknn或0,0qqkMqq00,1k0qtx)(0qMA n例2.3 设 为MA(1)序列：由它产生长度为n=300样本值 ，计算出前17个样本自相关函数为 ，计算出：tx) 1, 0(,35. 01NXtttt30021,xxx1710, ,rrr%3 .6875.68165160652. 021300121 2. 2. AICAIC准则定阶方法准

33、则定阶方法n给出模型阶数q的上界 ，对于 按前述的方法逐个拟合MA(m)模型。并给出白噪声方差 的估计量 ，定义AIC函数其中，n是样本个数，AIC(m)的最小值点(如不唯一，应取小的)称为MA(q)模型的AIC定阶。00Qm 22m0Q02, 1 , 0,2)(log)(QmnmmmAIC n例2.3的定阶问题，使用AIC准则，有q123456AIC2.8852.8912.8982.9022.8952.901 三三. . 拟合模型的检验拟合模型的检验n如果一段时间序列数据 符合(2.1)式，则当给定初始值 ，由(2.2)式计算出 ，它应当是独立序列 的一段样本值。故检验问题就转化为检验 是否

34、为独立序列的一段样本值的问题n方法： 检验和正态检验。nxxx,21q110,n,1tn,212 四四. . 建模例题：建模例题： 产生模型 的n=300个样本数据，建立模型。(1) 求出样本均值、样本自协方差函数、样本自相关函数),1 , 0(,5 . 01NXtttt00044071. 030011300nttxx n样本自相关函数-0.6-0.4-0.20.00.25101520253035RMX3 (2) 观察样本自相关函数为1步结尾，或使用前述的两种定阶方法，初步判定MA(1)(3) 使用第二小节的矩估计的解析方法可得(4) 检验：给出我们使用 检验，给出 ，计算出0389. 1)4

35、11 (2,6239. 0411221022111r0,6239. 001tttx2151 m)(300)(222212mm n 取值图 05101520253012345678910 11 12 13 14 15KATKA2 第三节第三节 ARMAARMA模型的拟合模型的拟合n根据数据序列 ，拟合以下ARMA(p,q)模型： (3.1)其中， 为独立同分布的序列，且 对一切st成立，参数和 满足平稳性和可逆性条件，且 与无公共根。nxxx,21qtqttptptttXXXX112211t22tE, 0tE,4tE),(21p),(1q)(u)(u, 0tsEx 一一. . 模型参数的估计模型

36、参数的估计1. 1. 矩估计方法矩估计方法步骤1. 的矩估计 ，满足如下方程： (3.2)其中， ，由(1.19)可知p元线性方程组。 qppqpqpqpqpqqqpqpqqqrrrrrrrrrrrr22112211211211qpkrrkk, 1 , 0, 记于是(3.2)可简写为若 满秩，则 (3.3)qpqpqpqqqpqqqqppqqqqprrrrrrrrrrrrb,212111,21, ,qpqpbqpqpb,1,)(qp, 步骤2. 和 满足以下的方程式 (3.4)式中其中， (3.4)式关于 的非线性代数方程组。当q=1,2可求出 显示解，当 ，可用数值解法。 2qkkyrqkq

37、kkqk1),(0),1 ()(1122212当当pjijikjikryr0,)(10q,1212,3q 2. 2. 近似极大似然估计方法近似极大似然估计方法n方法:取初始值 对于任意给定的一组参数 ,由(3.1)迭代算出相应值，即 (3.5)定义关于 的函数则，近似似然函数为),(n,21ptptttqtqttxxxx221111),(nkkS12),(),(21log22log),(222Snnl, 010pqxxx00pq 使得上式取到极大值的 ，称它们为 的近似极大似然估计，也称最小平方和估计。 当q=0，上述极值问题简化为Yule-Walker估计。 当p=0，上述极值问题与第三节的

38、近似极大似然估计方法一致。(3.1)式中的 的估计为), (),(2nS), (2 3. 3. 自回归逼近方法自回归逼近方法n基本思路：(1) 为数据 建立AR模型，取自回归阶数的上界 ，采用AIC定阶方法得到AR模型的阶数估计P和自回归系数 的估计. 。 (2) 计算残差 写出近似的ARMA(p,q)模型 nxxx,210nP ),(21pnPtxxxxptptttt, 1,2211qtqttptptttXXXX112211 (3) 对目标函数 (2.6)极小化,得到最小二乘估计 , 的最小二乘估计由下式定义nqPtPjqjjtjjtjtxxQ1112)(), (),(121qp2),(11

39、212qpQqPn n具体算法定义：pnnnpqPqPqPpqPqPqPqnnnPqPqPPqPqPnqPqPxxxxxxxxxExxx21211121211121,Xx 则目标函数(2.6)可写成，可解出最小二乘估计为相应地， 的估计为22|),(|)(XxXx,EEQxxXXXXxXXXTTTTTTTTTTEEEEEEEEE),(12) (12,QqPn 二二. . 模型阶数的估计模型阶数的估计1. 相关分析法用于ARMA模型的定阶n方法：(1) 给定初值 (一般取初值为零)将 的估计代入(3.2)递推得到残差估计 ,10pqxxxpq,0,n,1 (2) 作假设检验 来自于白噪声序列长度

40、为n的样本。 不是白噪声序列的长度为n的样本。令检验 等价于检验 是否来自于N(0,1)总体的k个独立抽样问题。nH,:10nH,:11kntkttknr11)(0Hn,1 a. 检验 的绝对值是否有68.3%左右小于1.b. 检验法：在 成立条件下，当n充分大时， 是k个相互独立N(0,1)随机变量，则 服从自由度为k的中心 分布，则以显著水平为 的否定域为kjnj, 2 , 1),(20Hkjnj, 2 , 1),(kjjknQ12)(22,kkQ 2. AIC准则方法n给定ARMA模型阶数的上界 和 。对于每一对(k,j)， ,计算AIC函数取 ，使 此时称 为ARMA模型的阶的估计，其

41、中 一般取 或 中的整数。0P0Q000 ,0QjPk,)(2),(log),(2njkjkjkAICqp , ),(min) , (0000qpAICqpAICQqPpqp , 00,QPn5/,10/nn n具有相合性的定阶准则BIC,使上式达最小的 为ARMA模型的阶。 中的整数。,log)(),(log),(2nnjkjkjkBICqp , ),(min) , (0000qpBICqpBICQqPp5/,10/,00nnQP 三三. . 拟合模型的检验拟合模型的检验nARMA模型的检验是检验其拟合残差序列是否为独立序列。n方法：取初值计算 的样本值，即 检验 是否为独立序列的样本值。0

42、010qpqxxxtntXxxxqtqtptptttt, 1,112211n,21 四四. . 例子：例子：kejian2kejian2由计算机产生模拟时间序列数据 。(1) 计算出样本均值、自相关函数、自协方差函数和偏相关函数) 1, 0(,6 . 03 . 011NYYttttt30021,yyy017102. 030011300nttyy n样本自相关函数-0.3-0.2-0.10.00.10.25101520253035RMY n样本偏相关函数-0.3-0.2-0.10.00.10.25101520253035PRMY (2) 取定阶数 由AIC准则：p=q=1(3) 估计参数：(4)

43、 检验 结论：数据符合ARMA(1,1)模型。99605. 0,75. 0,34375. 0211%47.7617/13057735. 03001 n -0.10-0.050.000.050.100.15246810121416第四节第四节 求和模型与季节模型求和模型与季节模型 的处理方法的处理方法一一. . 求和模型求和模型ARIMAARIMA的识别与拟合的识别与拟合1. 求和模型的识别方法n方法一 直接观察数据图形的方法：根据数据 画出数据曲线图，通过观察曲 线的形状，可初步判别是否需要拟合求和模型。nxxx,21 n例：数据1-3-2-1012350100150200250300Y n数

44、据2 -1001020309091929394959697989900TEMPX n数据3 0100200300400500600700495051525354555657585960X n数据40.0E+005.0E+071.0E+081.5E+082.0E+082.5E+082468101214161820X n方法二. 数据样本自相关函数分析法：当序列含有趋势项时，序列的样本自相关函数 的尾部不衰减到零值，特别地，所含趋势项为多项式时， 将近似于常数为1的序列。)(xk)(xk n例：序列105010015020025030050100 150 200 250 300 350 400 4

45、50 500X) 1 , 0(,5 . 0,2 . 021Nwwwtxtttttt n样本自相关函数：0.700.750.800.850.900.951.005101520253035404550RX n序列20200040006000800010000120001400050100 150 200 250 300 350 400 450 500Z) 1 , 0(,5 . 0,05. 01 . 0212Nwwwttxtttttt n序列2样本自相关函数0.91.05101520253035404550RZ n序列3 -4-20246850100 150 200 250 300

46、 350 400 450 500W) 1 , 0(,5 . 0,105. 0Nwwwexttttttt n序列3的样本自相关函数0.85101520253035404550RW n另外，还可从数据的来源判断使用求和模型的合理性。2. 2. 判断求和模型判断求和模型ARIMA(p,d,q)ARIMA(p,d,q)的阶数的阶数d dn对ARIMA(p,d,q)模型的研究焦点是对差分阶数d的判别。nd的判别方法：(1) 用动态数据 的实际背景来确定。若数据围绕着某条曲线变化，而此曲线是近似线性的，则判断差分阶数d=1,若此曲线可由二次多项式近似，则判断阶数d=2,

47、一般地，若该曲线可由d次t的多项式逼近，则可对原序列 作d次差分 ，而 可按平稳序列建模。nxxx,21txttdyx ty (2) 采用数据处理的方法：对原动态数据 分别作j次差分， ，连同原数据共有D+1套 动态数据，然后对每套数据求出样本自相关函数和样本 偏相关函数为 ,综合分析它 们的截尾性或拖尾性，最后判定为何种模型，再建立相 应的模型。 nxxx,21dj, 2 , 1djajjkkk, 2 , 1, 2. 2. 求和模型的拟合求和模型的拟合n步骤1：判断p值，对原数据进行d次差分运算，即 (4.1) 例：当d=2时， 为 的二次差分序列，即，nddtxBwtdt, 2, 1,)1

48、 (twtx2122)1 (tttttxxxxBw n步骤2：根据差分后的数据序列 按照前几节的方法，拟合AR、MA、ARMA模型，包括模型参数的估计以及对阶数p,q的估计，即 (4.2) 结合(4.1)和(4.2)，得到ARIMA(p,d,q)的拟合模型为， (4.3)n步骤3：对拟合求和模型 的检验：即是检验 是否符合(4.2)的模型，亦是对拟合后的残差进行白噪声检验.ndww,1ttBwB)()( ttdBxBB)()( )1 () , (qdpARIMAndww,1 n步骤3：对拟合求和模型 的检验：即是检验 是否符合(4.2)的模型，亦是 对拟合后的残差进行白噪声检验。ndww,1)

49、 , (qdpARIMA n例：某国1960年至1993年GNP平减指数的季度时间序列。要求对序列进行模型识别。(sample12)50100150200250657075808590X n )/log(1tttxxy-0.010.000.010.020.03657075808590Y n第一步：判断差分阶数d=1,对数据 进行一阶差分-0.015-0.010-0.0050.0000.0050.0100.0150.020657075808590Zty n第二步：对差分后序列进行ARMA模型拟合。观察样本自相关函数和偏相关函数，初步判断为AR模型。使用AIC定阶准则和最小二乘估计方法。判断阶数p

50、=2, ,即336226. 0,569612. 021003974. 0,1336226. 01569612. 0121ttttzzz n拟合后的残差图-0.02-0.010.000.010.02-0.02-0.010.000.010.02657075808590ResidualActualFitted n第三步：拟合模型的检验：采用正态检验， %3 .68%73.7211/8,087039. 0132/1-0.20-0.15-0.10-0.050.000.050.1060:160:361:161:362:162:3 n于是，拟合模型为：003974. 0)/log()336226. 0569

51、612. 0)(1 (121tttttttxxyyyyB 二二. . 季节模型的识别与拟合季节模型的识别与拟合n季节ARMA模型： 其中T是周期， 是某个ARMA(p,q)模型的特征多项式，实际问题中T经常的取值是4，7或12。tTtTwBXB)()()(),(uu 星期一二三四五六日一周二周n+1周1x2x7x8x9x14x3x4x5x6x10 x11x12x13xnx71nx72nx73nx74nx75nx76nx77 上述表中的每一列都可以看成一个时间序列，将数据(4.5)的第j列零均值化 (4.6) 其中 (4.7)首先，用数据(4.6)建立模型： (4.8)其中ttwBYB)()(7

52、7714271771427171)(1)(qqppBBBBBBBBntxyjtjtj, 1 , 0,77nttjjxn07) 1/(1 在相隔T步上为白噪声序列，而相隔小于T步时是相关的，即其次， 仍为平稳序列，故需对 建立ARMA(p,q)模型， (4.9)其中 对季节内外为白噪声序列，将(4.9)代入(4.8),得到季节ARMA 模型， (4.10)twTjkmTjkwEwEwjkt| , 0| , 0, 0twtwttBwB)()(t7),(),(QPqpttBBXBB)()()()(77 n季节模型(4.10)实际上是一个ARMA(p+7P，q+7Q)模型，只是其中又很多的系数是零。季

53、节模型的拟合方法：第一步：设 是数据 的样本自协方差函数，利用 拟合一个 模型： 要求这个模型通过模型检验；第二步：利用 拟合一个 模型：要求这个模型通过模型检验。于是， kr ty), 2 , 1 , 0(krik),(QPARMAttwBYB)()(771710, ,rrr) , (qpARMAttwBYB)()(77ttBBXBB)()()()( 77 n为了得到更精确的估计，可将模型看作疏系数的ARMA模型，使用前几节的ARMA模型参数的极大似然估计方法或最小二乘估计法估计模型(4.10)中的参数。 n 例：北京市1990.1-2000.12气温数据（sample6）-10010203

54、09091929394959697989900TEMPX n差分运算tempxBtempx)1 (112-6-4-20249091929394959697989900TEMPX1 n观察序列tempx1的样本自相关函数和偏相关函数，建模ls tempx ar(1) ar(2) ma(1) ma(2) sar(12) sma(12) -4-2024-100102030929394959697989900ResidualActualFitted n残差检验183340. 0119/2%,5 .9510/10-0.2-0.10.00.10.290:0190:0390:0590:0790:09 三三.

55、 . 乘积模型的拟合乘积模型的拟合n如果时间序列既具有趋势项又具有周期项，需采用乘积模型来拟合。在上例中如果每一列的数据需要经过差分后才能进行季节ARMA模型的拟合，模型将改写为， 称之为乘积模型 .ndDtyBBztDdt, 17,)1 ()1 (7ttDdBBXBBBB)()()1 ()1)()( 777),(),(QDPqdpARMA n实际问题中，d和D的取值一般很小，例如：D=0或1。季节模型 实际上是一个乘积模型 。n 一种简单的乘积模型： (4.11) 其中T为某一正整数，表示周期。 为某一平稳可逆的ARMA(p,q)序列。(1)(1)模型的拟合模型的拟合 在T和D已知时，首先对

56、 进行差分变换 (4.12),(),(QPqpARMA), 0 ,(), 0 ,(QPqpARMAttDTdWXBB)1 ()1 (tWtxnTDTDtXBytDTt, 2, 1,)1 (L 其中， 满足故，只需对 拟合形如(4.1)的求和模型，就可得到模型(4.11)的参数估计。(2) T2) T和和D D的取值判断的取值判断 a. T表示周期，它有较明显的物理背景，可根据数据的实际背景确定T的大小 b. D的确定，可使用逐步尝试的方法。即对D=1,2逐一尝试，并拟合(4.11),若模型检验通过，则确定该值为D的取值。tyttdwyB)1 (ty n例：航空客流量数据。01002003004

57、00500600700495051525354555657585960X n数据的预处理： (1) xx=log(x) (2) cx=xx-5.542193-1.0-0.50.00.51.0495051525354555657585960CX n模型的建立：对数据cx: )12()1 (12cxtcxcxBdcx-0.10.00.4495051525354555657585960DCX n对数据dcx进行一次差分运算-0.15-0.10-0.050.000.050.100.15495051525354555657585960DDCX) 1(dcxdcxddcx n对数据ddc

58、x建立模型：使用最小二乘估计方法，得到13312211tttttbbbddcx035247. 0,468081. 0,747640. 0,373970. 0321bbbls ddcx ma(1) ma(2) ma(3) n模型的残差检验：-0.20-0.15-0.10-0.050.000.050.100.1549:0149:0349:0549:0749:0949:11%4 .68%73.7211/8 n结论：客流量的模型为035247. 0468081. 0747640. 0373970. 0)542193. 5)(log(1312112tttttx n另一种方法建模：观察序列ddcx的样本自

59、相关函数和偏相关函数。ls d(log(x),1,12) ar(1) ar(2) ar(3) ma(1) sar(12) sma(12)ls ddcx ar(1) ar(2) ar(3) ma(1) sar(12) sma(12)ls d(log(x),1,12) ar(1) ar(3) ma(1) sar(12) sma(12)ls d(log(x),1,12) ar(3) ma(1) sar(12) sma(12)-0.10-0.050.000.050.100.1549:0149:0349:0549:0749:0949:11ls d(log(x),1,12) ar(3) ma(1) sar

60、(12) sma(12) 1. 设 为零均值平稳序列，由它的长度为N=100的样本算得样本自相关函数 及样本偏相关函数 的前6个数值如下又知 ，试求：(1) 为哪种模型，并说明理由。(2) 对模型参数和白噪声方差给出矩估计。(3) 判断所建立的模型是否具有平稳性（或可逆性），并给出模型的传递形式（或逆转形式）txkkka34. 30tx k123456-0.8000.670-0.1580.390-0.3100.221-0.8000.0850.112-0.046-0.0610.038kkka 2. 全国城镇居民储蓄额年数据序列 的建模。下表给出1952年至1991年储蓄额年数据(亿元)tx8.6