《算法大全第24章 时间序列模型.pdf》由会员分享，可在线阅读，更多相关《算法大全第24章 时间序列模型.pdf（29页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、-280-第二十四章第二十四章 时间序列模型时间序列模型 时间序列是按时间顺序排列的、随时间变化且相互关联的数据序列。分析时间序列的方法构成数据分析的一个重要领域，即时间序列分析。时间序列根据所研究的依据不同，可有不同的分类。1按所研究的对象的多少分，有一元时间序列和多元时间序列。2按时间的连续性可将时间序列分为离散时间序列和连续时间序列两种。3按序列的统计特性分，有平稳时间序列和非平稳时间序列。如果一个时间序列的概率分布与时间t无关，则称该序列为严格的（狭义的）平稳时间序列。如果序列的一、二阶矩存在，而且对任意时刻t满足：（1）均值为常数（2）协方差为时间间隔的函数。则称该序列为宽平稳时间序
2、列，也叫广义平稳时间序列。我们以后所研究的时间序列主要是宽平稳时间序列。4按时间序列的分布规律来分，有高斯型时间序列和非高斯型时间序列。1 确定性时间序列分析方法概述 时间序列预测技术就是通过对预测目标自身时间序列的处理，来研究其变化趋势的。一个时间序列往往是以下几类变化形式的叠加或耦合。（1）长期趋势变动。它是指时间序列朝着一定的方向持续上升或下降，或停留在某一水平上的倾向，它反映了客观事物的主要变化趋势。（2）季节变动。（3）循环变动。通常是指周期为一年以上，由非季节因素引起的涨落起伏波形相似的波动。（4）不规则变动。通常它分为突然变动和随机变动。通常用tT表示长期趋势项，tS表示季节变动
3、趋势项，tC表示循环变动趋势项，tR表示随机干扰项。常见的确定性时间序列模型有以下几种类型：（1）加法模型 tttttRCSTy+=（2）乘法模型 tttttRCSTy=（3）混合模型 ttttRSTy+=tttttRCTSy+=其中ty是观测目标的观测记录，0)(=tRE，22)(=tRE。如果在预测时间范围以内，无突然变动且随机变动的方差2较小，并且有理由认为过去和现在的演变趋势将继续发展到未来时，可用一些经验方法进行预测。2 移动平均法 移动平均法是根据时间序列资料逐渐推移，依次计算包含一定项数的时序平均数，以反映长期趋势的方法。当时间序列的数值由于受周期变动和不规则变动的影响，起伏较大
4、，不易显示出发展趋势时，可用移动平均法，消除这些因素的影响，分析、预测序-281-列的长期趋势。移动平均法有简单移动平均法，加权移动平均法，趋势移动平均法等。2.1 简单移动平均法 设观测序列为Tyy,1L，取移动平均的项数TN。一次简单移动平均值计算公式为：)(111)1(+=NttttyyyNML )(1)(1)(1)1(11NtttNttNttyyNMyyNyyN+=+=L （1）当预测目标的基本趋势是在某一水平上下波动时，可用一次简单移动平均方法建立预测模型：)(11)1(1+=NttttyyNMyL，L,1,+=NNt，（2）其预测标准误差为：NTyySTNttt=+=12)(，（3
5、）最近N期序列值的平均值作为未来各期的预测结果。一般N取值范围：2005 N。当历史序列的基本趋势变化不大且序列中随机变动成分较多时，N的取值应较大一些。否则N的取值应小一些。在有确定的季节变动周期的资料中，移动平均的项数应取周期长度。选择最佳N值的一个有效方法是，比较若干模型的预测误差。预测标准误差最小者为好。例 1 某企业 1 月11 月份的销售收入时间序列如表 1 示。试用一次简单滑动平均法预测第 12 月份的销售收入。表 1 企业销售收入 月份t 1 2 3 4 5 6 销售收入ty 533.8 574.6 606.9 649.8 705.1 772.0 月份t 7 8 9 10 11
6、 销售收入ty 816.4 892.7 963.9 1015.1 1102.7 解：分别取5,4=NN的预测公式 4321)1(1+=tttttyyyyy，11,5,4L=t 54321)2(1+=ttttttyyyyyy，11,5L=t 当4=N时，预测值993.6)1(12=y，预测的标准误差为 150.5411)(1152)1(1=tttyyS 当5=N时，预测值182.4)2(12=y，预测的标准误差为 -282-958.2511)(1162)2(2=tttyyS 计算结果表明，4=N时，预测的标准误差较小，所以选取4=N。预测第 12 月份的销售收入为 993.6。计算的 Matla
7、b 程序如下：clc,clear y=533.8 574.6 606.9 649.8 705.1 772.0 816.4 892.7 963.9 1015.1 1102.7;m=length(y);n=4,5;%n 为移动平均的项数 for i=1:length(n)%由于 n 的取值不同，yhat 的长度不一致，下面使用了细胞数组 for j=1:m-n(i)+1 yhati(j)=sum(y(j:j+n(i)-1)/n(i);end y12(i)=yhati(end);s(i)=sqrt(mean(y(n(i)+1:m)-yhati(1:end-1).2);end y12,s 简单移动平均
8、法只适合做近期预测，而且是预测目标的发展趋势变化不大的情况。如果目标的发展趋势存在其它的变化，采用简单移动平均法就会产生较大的预测偏差和滞后。2.2 加权移动平均法 在简单移动平均公式中，每期数据在求平均时的作用是等同的。但是，每期数据所包含的信息量不一样，近期数据包含着更多关于未来情况的信心。因此，把各期数据等同看待是不尽合理的，应考虑各期数据的重要性，对近期数据给予较大的权重，这就是加权移动平均法的基本思想。设时间序列为LL,21tyyy；加权移动平均公式为 NNtNttwwwwywywywM+=+LL211221，Nt （4）式中twM为t期加权移动平均数；iw为1+ity的权数，它体现
9、了相应的ty在加权平均数中的重要性。利用加权移动平均数来做预测，其预测公式为 twtMy=+1 （5）即以第t期加权移动平均数作为第1+t期的预测值。例 2 我国 19791988 年原煤产量如表 2 所示，试用加权移动平均法预测 1989 年的产量。表 2 我国原煤产量统计数据及加权移动平均预测值表 年份 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 原煤产量ty 6.35 6.20 6.22 6.66 7.15 7.89 8.72 8.94 9.28 9.8 三年加权移动平均预测值 6.2356.43676.83177.43838.18
10、17 8.6917 9.0733 -283-相对误差（）6.38 9.98 13.4114.7 8.48 6.34 7.41 解 取1,2,3321=www，按预测公式 12323211+=+ttttyyyy 计算三年加权移动平均预测值，其结果列于表 2 中。1989 年我国原煤产量的预测值为（亿吨）48.9694.828.928.931989=+=y 这个预测值偏低，可以修正。其方法是：先计算各年预测值与实际值的相对误差，例如1982 年为%38.666.6235.666.6=将相对误差列于表 2 中，再计算总的平均相对误差。%5.9%100)44.5889.521(%1001=ttyy 由
11、于总预测值的平均值比实际值低%5.9，所以可将 1989 年的预测值修正为 4788.10%5.9148.9=计算的 MATLAB 程序如下：y=6.35 6.20 6.22 6.66 7.15 7.89 8.72 8.94 9.28 9.8;w=1/6;2/6;3/6;m=length(y);n=3;for i=1:m-n+1 yhat(i)=y(i:i+n-1)*w;end yhat err=abs(y(n+1:m)-yhat(1:end-1)./y(n+1:m)T_err=1-sum(yhat(1:end-1)/sum(y(n+1:m)y1989=yhat(end)/(1-T_err)在
12、加权移动平均法中，tw的选择，同样具有一定的经验性。一般的原则是：近期数据的权数大，远期数据的权数小。至于大到什么程度和小到什么程度，则需要按照预测者对序列的了解和分析来确定。2.3 趋势移动平均法 简单移动平均法和加权移动平均法，在时间序列没有明显的趋势变动时，能够准确反映实际情况。但当时间序列出现直线增加或减少的变动趋势时，用简单移动平均法和加权移动平均法来预测就会出现滞后偏差。因此，需要进行修正，修正的方法是作二次移动平均，利用移动平均滞后偏差的规律来建立直线趋势的预测模型。这就是趋势移动平均法。一次移动的平均数为 -284-)(111)1(+=NttttyyyNML 在一次移动平均的基
13、础上再进行一次移动平均就是二次移动平均，其计算公式为 )(1)(1)1()1()2(1)1(1)1()2(NtttNtttMMNMMMNM+=+=L （6）下面讨论如何利用移动平均的滞后偏差建立直线趋势预测模型。设时间序列ty从某时期开始具有直线趋势，且认为未来时期也按此直线趋势变化，则可设此直线趋势预测模型为 TbayttTt+=+，L,2,1=T （7）其中t为当前时期数；T为由t至预测期的时期数；ta为截距；tb为斜率。两者又称为平滑系数。现在，我们根据移动平均值来确定平滑系数。由模型（7）可知 ttya=tttbyy=1 tttbyy22=ttNtbNyy)1(1=+所以 tttttt
14、tttNttttbNyNbNNyNbNybyyNyyyM21)1(21)1()(11)1(=+=+=+=+LLL 因此 tttbNMy21)1(=（8）由式（7），类似式（8）的推导，可得 tttbNMy21)1(11=（9）所以 tttttbMMyy=)1(1)1(1 （10）类似式（8）的推导，可得 tttbNMM21)2()1(=（11）于是，由式（8）和式（11）可得平滑系数的计算公式为 =)(122)2()1()2()1(ttttttMMNbMMa （12）例 3 我国 19651985 年的发电总量如表 3 所示，试预测 1986 年和 1987 年的发电总量。-285-表 3 我
15、国发电量及一、二次移动平均值计算表 年份 t 发电总量 yt 一次移动平均，N6 二次移动平均，N61965 1 676 1966 2 825 1967 3 774 1968 4 716 1969 5 940 1970 6 1159 848.3 1971 7 1384 966.3 1972 8 1524 1082.8 1973 9 1668 1231.8 1974 10 1688 1393.8 1975 11 1958 1563.5 1181.1 1976 12 2031 1708.8 1324.5 1977 13 2234 1850.5 1471.9 1978 14 2566 2024.2
16、1628.8 1979 15 2820 2216.2 1792.8 1980 16 3006 2435.8 1966.5 1981 17 3093 2625 2143.4 1982 18 3277 2832.7 2330.7 1983 19 3514 3046 2530 1984 20 3770 3246.7 2733.7 1985 21 4107 3461.2 2941.2 解 由散点图 1 可以看出，发电总量基本呈直线上升趋势，可用趋势移动平均法来预测。051015202550010001500200025003000350040004500 图 1 原始数据散点图 取6=N，分别计算一次
17、和二次移动平均值并列于表 3 中。2.3461)1(21=M，2.2941)2(21=M 再由公式（12），得 3981.12)2(21)1(2121=MMa 208)(162)2(21)1(2121=MMb 于是，得21=t时直线趋势预测模型为 TyT2081.398121+=+预测 1986 年和 1987 年的发电总量为 1.4192121221986=+yyy 1.4397221231987=+yyy 计算的 MATLAB 程序如下：-286-clc,clear load y.txt%把原始数据保存在纯文本文件 y.txt 中 m1=length(y);n=6;%n 为移动平均的项数
18、for i=1:m1-n+1 yhat1(i)=sum(y(i:i+n-1)/n;end yhat1 m2=length(yhat1);for i=1:m2-n+1 yhat2(i)=sum(yhat1(i:i+n-1)/n;end yhat2 plot(1:21,y,*)a21=2*yhat1(end)-yhat2(end)b21=2*(yhat1(end)-yhat2(end)/(n-1)y1986=a21+b21 y1987=a21+2*b21 趋势移动平均法对于同时存在直线趋势与周期波动的序列，是一种既能反映趋势变化，又可以有效地分离出来周期变动的方法。3 指数平滑法 一次移动平均实际
19、上认为最近N期数据对未来值影响相同，都加权N1；而N期以前的数据对未来值没有影响，加权为 0。但是，二次及更高次移动平均数的权数却不是N1，且次数越高，权数的结构越复杂，但永远保持对称的权数，即两端项权数小，中间项权数大，不符合一般系统的动态性。一般说来历史数据对未来值的影响是随时间间隔的增长而递减的。所以，更切合实际的方法应是对各期观测值依时间顺序进行加权平均作为预测值。指数平滑法可满足这一要求，而且具有简单的递推形式。指数平滑法根据平滑次数的不同，又分为一次指数平滑法、二次指数平滑法和三次指数平滑法等，分别介绍如下。3.1 一次指数平滑法 1预测模型 设时间序列为LL,21tyyy，为加权
20、系数，100.00001 Terr=;for j=N+1:m-1 yhat(j)=w*yt(j-1:-1:j-N);err=yt(j)-yhat(j);Terr=Terr,abs(err);w=w+2*k*err*yt(j-1:-1:j-N);end Terr=max(Terr);end w,yhat 5.2 kN,值和初始权数的确定 在开始调整权数时，首先要确定权数个数N和学习常数k。一般说来，当时间序列的观测值呈季节变动时，N应取季节性长度值。如序列以一年为周期进行季节变动时，若数据是月度的，则取12=N，若季节是季度的，则取4=N。如果时间序列无明显的周期变动，则可用自相关系数法来确定，
21、即取N为最高自相关系数的滞后时期。k的取值一般可定为N/1，也可以用不同的k值来进行计算，以确定一个能使S最小的k值。初始权数的确定也很重要，如无其它依据，也可用N/1作为初始权系数用，即 ),2,1(1NiNwiL=自适应滤波法有两个明显的优点：一是技术比较简单，可根据预测意图来选择权数的个数和学习常数，以控制预测。也可以由计算机自动选定。二是它使用了全部历史数据来寻求最佳权系数，并随数据轨迹的变化而不断更新权数，从而不断改进预测。由于自适应滤波法的预测模型简单，又可以在计算机上对数据进行处理，所以这种预测方法应用较为广泛。6 趋势外推预测方法 趋势外推法是根据事物的历史和现时资料，寻求事物
22、发展规律，从而推测出事物未来状况的一种比较常用的预测方法。利用趋势外推法进行预测，主要包括六个阶段：（a）选择应预测的参数；（b）收集必要的数据；（c）利用数据拟合曲线；（d）趋势外推；（e）预测说明；（f）研究预测结果在进行决策中应用的可能性。趋势外推法常用的典型数学模型有：指数曲线、修正指数曲线、生长曲线、包络曲线等。6.1 指数曲线法 一般来说，技术的进步和生产的增长，在其未达饱和之前的新生时期是遵循指数曲线增长规律的，因此可以用指数曲线对发展中的事物进行预测。指数曲线的数学模型为 Kteyy0=（35）其中系数0y和K值由历史数据利用回归方法求得。对式（35）取对数可得 Ktyy+=0
23、lnln （36）-297-令 yYln=，0lnyA=则 KtAY+=其中KA,可以用最小二乘法求得。6.2 修正指数曲线法 利用指数曲线外推来进行预测时，存在着预测值随着时间的推移会无限增大的情况。这是不符合客观规律的。因为任何事物的发展都是有一定限度的。例如某种畅销产品，在其占有市场的初期是呈指数曲线增长的，但随着产品销售量的增加，产品总量接近于社会饱和量时。这时的预测模型应改用修正指数曲线。ttabKy+=（37）在此数学模型中有三个参数aK,和b要用历史数据来确定。修正指数曲线用于描述这样一类现象。（1）初期增长迅速，随后增长率逐渐降低。（2）当0K，0a，10 b时，t，0tab，
24、即Kyt。当K值可预先确定时，采用最小二乘法确定模型中的参数。而当K值不能预先确定时，应采用三和法。把时间序列的n个观察值等分为三部分，每部分有m期，即mn3=。第一部分：myyy,21L；第二部分：mmmyyy221,L+；第三部分：mmmyyy32212,L+令每部分的趋势值之和等于相应的观察值之和，由此给出参数估计值。三和法步骤如下：记观察值的各部分之和+=+=mmttmmttmttySySyS312321211,（38）且 +=+=+=+=+=+=+=+=+=+=)1()()1()()1()(121231231231212121212111mmmmttmmttmmmmttmmttmmt
25、tmttbbbabmKabKySbbbabmKabKySbbbabmKabKySLLL （39）由于 1)1)(1(12=+mmbbbbbL （40）则根据(11.35)式，得 -298-+=+=+=+111112311211bbabmKSbbabmKSbbabmKSmmmmm （41）由（41）式，解得=1)1(1)1(1)(121211223bbabSmKbbbSSaSSSSbmmm （42）至此三个参数全部确定了，于是就可以用式（37）进行预测。值得注意的是，并不是任何一组数据都可以用修正指数曲线拟合。采用前应对数据进行检验，检验方法是看给定数据的逐期增长量的比率是否接近某一常数b。即
26、byyyytttt+11 （43）例 8 根据统计资料，某厂收音机连续 15 年的销售量如表 10。表 11 某厂收音机销售量 时间(年)1969 1970 1971 1972 1973 1974 1975 1976 销售量(万部)42.1 47.5 52.7 57.7 62.5 67.1 71.5 75.7 时间(年)1977 1978 1979 1980 1981 1982 1983 销售量(万部)79.8 83.7 87.5 91.1 94.6 97.9 101.1 试用修正指数曲线预测 1986 年的销售量。解 经计算可知 0.9762 0.9429,11+ttttyyyy 可以认定这
27、组数据可以采用修正指数曲线拟合。现将以上 15 个数据分为三部分，每部分 5 个数据，即15=n，5=m，并以 1969 年作为开始年份1=t。根据式（38），得 5.2621=S，8.3772=S，2.4723=S 再由（42）式，得 9608.0=b，2063.143=a，7162.179=K 故修正指数曲线的数学模型为 ty9608.02063.1437162.179=（44）预测 1986 年的产量时，18119691986=+=t。所以 -299-1109608.02063.1437162.179181986=y（万部）计算的 MATLAB 程序如下：function chanlia
28、ng clc,clear global a b k load xsh.txt%原始数据存放在纯文本文件 xsh.txt 中 yt=xsh;n=length(yt);m=n/3 cf=diff(yt);for i=1:n-2 bzh(i)=cf(i+1)/cf(i);end range=minmax(bzh)s1=sum(yt(1:m),s2=sum(yt(m+1:2*m),s3=sum(yt(2*m+1:end)b=(s3-s2)/(s2-s1)(1/m)a=(s2-s1)*(b-1)/(b*(bm-1)2)k=(s1-a*b*(bm-1)/(b-1)/m y=yuce(1:18)%*%定义预
29、测函数%*function y=yuce(t)global a b k y=k+a*b.t;6.3 Compertz 曲线 曲线的一般形式 tbtKay=，10,10,0baK （45）采用 Compertz 曲线前应对数据进行检验，检验方法是看给定数据的对数逐期增长量的比率是否接近某一常数b。即 byyyytttt+11lnlnlnln （46）Compertz 曲线用于描述这样一类现象：初期增长缓慢，以后逐渐加快。当达到一定程度后，增长率又逐渐下降。参数估计方法如下：式（46）两边取对数，得 ttbaKy)(logloglog+=（47）记 aaKKyyttln,ln,ln=得 ttbaK
30、y+=仿照修正指数曲线的三和法估计参数，令+=+=mmttmmttmttySySyS312321211,（48）其中ttyyln=。则类似式（42），得 -300-=1)1(1)1(1)(121211223bbbaSmKbbbSSaSSSSbmmm （49）例 9（续例 8）根据表 10 的数据，试确定收音机销售量的 Gompertz 曲线方程，求出各年收音机销售量的趋势值，并预测 1986 年的销售量。解：已知15=n，5=m，根据式（48），得 19.75581=S，21.60942=S，22.73333=S 再由式（49），得 0.9048=b 1.2588=a，0.284=a 4.89
31、29=K，133.3341=K 从而收音机销售量的 Compertz 曲线方程为 tty9048.0284.03341.133=将18=t代入方程，得 1986 年收音机销售量的预测值为 108.31431986=y（万部）计算的 MATLAB 程序如下：function chanliang2 clc,clear global a b k load xsh.txt%原始数据存放在纯文本文件 xsh.txt 中 yt=log(xsh);n=length(yt);m=n/3;s1=sum(yt(1:m),s2=sum(yt(m+1:2*m),s3=sum(yt(2*m+1:end)b=(s3-s2
32、)/(s2-s1)(1/m)a=(s2-s1)*(b-1)/(b*(bm-1)2)k=(s1-a*b*(bm-1)/(b-1)/m a=exp(a)k=exp(k)y=yuce(1:18)%*%定义预测函数%*function y=yuce(t);global a b k y=k*a.(b.t);6.4 Logistic 曲线（生长曲线）生物的生长过程经历发生、发展到成熟三个阶段，在三个阶段生物的生长速度是不一样的，例如南瓜的重量增长速度，在第一阶段增长的较慢，在发展时期则突然加快，而到了成熟期又趋减慢，形成一条 S 形曲线，这就是有名的 Logistic 曲线（生长曲线），-301-很多事物
33、，如技术和产品发展进程都有类似的发展过程，因此 Logistic 曲线在预测中有相当广泛的应用。Logistic 曲线的一般数学模型是 )1(Lyrydtdy=（50）式中y为预测值，L为y的极限值，r为增长率常数，0r。解此微分方程得 rtceLy+=1 （51）式中c为常数。下面我们记 Logistic 曲线的一般形式为 ttabKy+=1，0K，0a，10 b （52）检验能否使用 Logistic 曲线的方法，是看给定数据倒数的逐期增长量的比率是否接近某一常数b。即 byyyytttt+11/1/1/1/1 （53）Logistic 曲线中参数估计方法如下：作变换 ttyy1=得 tt
34、abKy+=仿照修正指数曲线的三和法估计参数，令 +=+=mmttmmttmttySySyS312321211,（54）则类似式（42），得=1)1(1)1(1)(121211223bbabSmKbbbSSaSSSSbmmm （55）例 10（续例 8）根据表 10 的数据，试确定收音机销售量的 Logistic 曲线方程，求出各年收音机销售量的趋势值，并预测 1986 年的销售量。解：已知15=n，5=m，根据式（54），得 0.09711=S，0.06662=S，0.05313=S 再由式（55），得 0.8493=b，0.0174=a，0.0085=K -302-从而收音机销售量的 Lo
35、gistic 曲线方程为 tty8493.00174.00085.01+=将18=t代入方程，得 1986 年收音机销售量的预测值为 106.39811986=y 计算的 MATLAB 程序如下：function chanliang3 clc,clear global a b k load xsh.txt%原始数据存放在纯文本文件 xsh.txt 中 yt=1./xsh;n=length(yt);m=n/3;s1=sum(yt(1:m),s2=sum(yt(m+1:2*m),s3=sum(yt(2*m+1:end)b=(s3-s2)/(s2-s1)(1/m)a=(s2-s1)*(b-1)/(b
36、*(bm-1)2)k=(s1-a*b*(bm-1)/(b-1)/m y=yuce(1:18)%*%定义预测函数%*function y=yuce(t);global a b k y=1./(k+a*b.t);6.5 趋势线的选择 趋势线的选择有以下几种方式。1由散点图选择趋势线。2由数据本身的取值规律选择趋势线。3比较预测标准误差大小 =niiiyynS12)(1 （56）当有几种趋势线可供选择时，应选择S最小的趋势线。7 平稳时间序列模型 这里的平稳是指宽平稳，其特性是序列的统计特性不随时间的平移而变化，即均值和协方差不随时间的平移而变化。下面自回归模型（Auto Regressive Mo
37、del）简称 AR 模型，移动平均模型（Moving Average Model）简称 MA 模型，自回归移动平均模型（Auto Regressive Moving Average Model）简称 ARMA 模型。下面的tX为零均值（即中心化处理的）平稳序列。（1）一般自回归模型 AR(n)假设时间序列tX仅与ntttXXX,21L有线性关系，而在ntttXXX,21L已 知 条 件 下，tX与),2,1(L+=nnjXjt无 关，ta是 一 个 独 立 于ntttXXX,21L的白噪声序列，),0(2atNa。tntntttaXXXX+=L2211 -303-上式还可以表示为 ntnttt
38、tXXXXa=L2211 可见，)(ARn系统的响应tX具有n阶动态性。)(ARn模型通过把tX中的依赖于ntttXXX,21L的部分消除掉之后，使得具有n阶动态性的序列tX转化为独立的序列ta。因此，拟合)(ARn模型的过程也就是使相关序列独立化的过程。（2）移动平均模型)(MA m)(ARn系 统 的 特 征 是 系 统 在t时 刻 的 响 应tX仅 与 其 以 前 时 刻 的 响 应ntttXXX,21L有关，而与其以前时刻进入系统的扰动无关。如果一个系统在t时刻的响应tX，与其以前时刻L,2,1 tt的响应L,21ttXX无关，而与其以前时刻mttt,2,1L进入系统的扰动mtttaa
39、a,21L存在着一定的相关关系，那么，这一类系统为)(MAm系统。mtmttttaaaaX=L2211（3）自回归移动平均模型 一个系统，如果它在时刻t的响应tX，不仅与其以前时刻的自身值有关，而且还与其以前时刻进入系统的扰动存在一定的依存关系，那么，这个系统就是自回归移动平均系统。),(ARMAmn模型为 mtmttntnttaaaXXX=LL1111 对于平稳系统来说，由于 AR、MA、),(ARMAmn模型都是)1,(ARMAnn模型的特例，我们以)1,(ARMAnn模型为一般形式来建立时序模型。8 ARMA 模型的特性 在时间序列的时域分析中，线性差分方程是极为有效的工具。事实上，任何
40、一个ARMA 模型都是一个线性差分方程。8.1 AR(1)系统的格林函数 格林函数就是描述系统记忆扰动程度的函数。AR(1)模型为 tttaXX=11 设)(1tyXt=，则有 tatyty=+)()1(1 （57）显然是一个一阶非齐次差分方程。由于 ttttttaaXaXX+=+=)(121111 tttaaX+=11221 L=+=ttttaaaX11221331 依次递推下去，可得到=01jjtjtaX （58）-304-显然（58）式是差分方程（57）的解。方程解的系数函数j1客观地描述了该系统的动态性，故这个系统函数就叫做记忆函数，也叫格林函数。若用jG表示，则 AR(1)模型的格林
41、函数可以表示为 jjG1=这样（58）式可等价地写成：=0jjtjtaGX （59）上式也可写成=tkkkttaGX 这里1010=G。定义后移算子B，L,221=ttttXXBXBX，这样，AR(1)可写成 ttaXB=)1(1 它的解为 tttaBBaBX)1(1122111L+=+=022111jjtjtttaGaaaL 由于格林函数就是差分方程解的系数函数，格林函数的意义可概括如下：（1）jG是前j个时间单位以前进入系统的扰动jta对系统现在行为（响应）影响的权数。（2）jG客观地刻画了系统动态响应衰减的快慢程度。（3）对于一个平稳系统来说，在某一时刻由于受到进入系统的扰动ta的作用，
42、离开其平衡位置（即平均数零），jG描述系统回到平衡位置的速度，1的值较小，速度较快；1的值较大，回复的速度就较慢。8.2 )1,2(ARMA系统的格林函数（1）)1,2(ARMA系统的格林函数的隐式)1,2(ARMA模型是一个二阶非齐次差分方程 112211=tttttaaXXX 设该二阶非齐次差分方程的解为=0jjtjtaGX （60）为方便起见，可用B算子：ttaBXBB)1()1(1221=（61）tjjjtaBGX=0 （62）-305-把（62）式代入（61），比较两边B的同次幂的系数得 10=G，111=G，L,4,3,2211=+=jGGGjjj（2）)1,2(ARMA系统的格林
43、函数的显式)1,2(ARMA模型实质上是一个二阶非齐次差分方程：112211=tttttaaXXX 欲求其解，必须先求出其相应的齐次差分方程的通解。齐次方程对应的特征方程为 0212=，特征根为 2422111+=，2422112+=齐次差分方程的通解为 jjjccG2211+=其中1c和2c是任意常数，其值由初始条件唯一地确定。这里的初始条件为：=11111101GGG 于是有 =+=+=11221112101ccGccG 而121=+，即+=+=+1212211211cccc 解之，得 21111=c，12122=c 则)1,2(ARMA系统的格林函数为：jjjG2121212111+=8
44、.3 逆函数和可逆性 前面的格林函数，把tX表示为过去ta对tX的影响，或者说系统对过去ta的记忆性，也就是用一个 MA 模型来逼近tX的行为。平稳序列tX的这种表达形式称为tX的“传递形式”。同样我们也可以用过去的tX的一个线性组合来逼近系统现在时刻的行为。即=+=+=12211jtjtjttttaXIaXIXIXL 我们把这种表达形式称为tX的“逆转形式”。其中的系数函数)1(0=IIj称为逆函-306-数。可见它是一个无穷阶的自回归模型。一个过程是否具有逆转形式，也就是说逆函数是否存在的性质，通常称为过程是否具有可逆性，如果一个过程可以用一个无限阶的自回归模型逼近，即逆函数存在，我们就称
45、该过程具有可逆性，否则，就是不可逆的。对于 AR（2）模型 ttttaXXX+=2211 有 11=I，22=I，L,4,3,0=jIj 可见，所谓可逆性，是指移动平均模型可以用 AR 模型表示。MA(1)模型：ttaBX)1(1=那么=+=+=1122111)1(1jjtjttttXXXBBBXaL 即=+=11)(jtjtjtaXX 可见，jjI1=，显然，只有1|1时，才有j，0jI，故 MA（1）的可逆性条件为 1|1=1jIIjj，于是 231230IIII=1。注意：以 AR(3)中的321,III替代 ARMA(2,1)中的321,III是一种近似代替。通-307-过这种方法求得
46、的1的绝对值若大于 1，则取其倒数作为初始值，以满足可逆性条件。知道了321,III及1，再用下式来确定 ARMA(2,1)模型中的21,：111+=I；21122II+=。（3）以（2）中得到的121,为初始值，利用非线性最小二乘法得到121,的终值及置信区间，并且求出残差平方和(RSS)。31+=nn，拟合 ARMA)12,2(nn模型 其基本步骤与 2 类似。4用 F 准则检验模型的适用性。若F检验显著，则转入第 2 步。若F检验不显著，转入第 5 步。对于 ARMA 模型的适用性检验的实际就是对ta的独立性检验。检验ta的独立性的一个简便而有效的办法是拟合更高阶的模型。若更高阶模型的残
47、差平方和有明显减少，就意味着现有模型的ta不是独立的，因而模型不适用；若更高阶模型的残差平方和没有明显减少，同时更高阶模型中的附加参数的值也很小（其置信区间包含 0），则可认为该模型是适用的。具体的检验准则如下。设 有 模 型),(ARMA11mn和),(ARMA22mn，1212,mmnn。假 设),(ARMA110mnA=模型的残差ta之平方和，),(ARMA221mnA=模型的残差ta之平方和，N是采集数据的数目，则检验准则为：),(001=NsFNAsAAF，其中22mn+=，)(1122mnmns+=。若这样得到的F值超过由F分布查表所得的在 5%置信水平上的),(NsF值，那么由)
48、,(ARMA11mn模型改变为),(ARMA22mn时，残差平方和的改善是显著的，因而拒绝关于模型),(ARMA11mn的适用性假设；F值低于查表所得之值，就可以认为在该置信水平上这个模型是适用的。5检查122,nn的值是否很小，其置信区间是否包含零。若不是，则适用的模型就是)12,2ARMA(nn。若122,nn很小，且其置信区间包含零，则拟合)22,12ARMA(nn。6利用 F 准则检验模型)12,2ARMA(nn和)22,12ARMA(nn，若 F 值不显著，转入第 7 步；若 F 值显著，转入第 8 步。7舍弃小的 MA 参数，拟合22 nm的模型),12ARMA(mn，并用 F 准
49、则进行检验。重复这一过程，直到得出具有最小参数的适用模型为止。8舍弃小的 MA 参数，拟合12 nm的模型),2ARMA(mn，并用 F 准则进行检验。重复这一过程，直到得出具有最小参数的适用模型为止。习题二十四 1我国 19741981 年布的产量如表 11 所示。表 11 19741981 年布的产量 年份 1974 1975 1976 1977 1978 1979 1980 1981 产量(亿米)80.8 94.0 88.4 101.5 110.3 121.5 134.7 142.7 -308-（1）试用趋势移动平均法（取3=N），建立布的年产量预测模型。（2）分别取3.0=，6.0=，
50、7.873321)2(0)1(0=+=yyySS，建立布的直线指数平滑预测模型。（3）计算模型拟合误差，比较 3 个模型的优劣。（4）用最优的模型预测 1982 年和 1985 年布的产量。219601982 年全国社会商品零售额如表 12 所示（单位：亿元）。表 12 全国社会商品零售额数据 年份 1960 1961 1962 1963 1964 1965 1966 1967 零售总额 696.9 607.7 604 604.5 638.2 670.3 732.8 770.5 年份 1968 1969 1970 1971 1972 1973 1974 1975 零售总额 737.3 801.