应用时间序列分析实验报告.docx
《应用时间序列分析实验报告.docx》由会员分享，可在线阅读，更多相关《应用时间序列分析实验报告.docx（27页珍藏版）》请在冰点文库上搜索。
应用时间序列分析实验报告
应用时间序列分析实验报告
学院名称 理学院
专业班级 应用统计学14-2
学生姓名 张艳雪
学号 201411081051
齐鲁工业大学实验报告成绩
课程名称《应用时间序列分析实验》指导教师黄玉林实验日期2017.6.30
院（系）理学院专业班级统计14-2实验地点机电楼C428
学生姓名张艳雪学号201411081051同组人无
实验项目名称ARIMA模型、确定性分析法,多元时间序列建模
一、实验目的和要求
1.熟悉非平稳序列的确定性分析法：
趋势分析、季节效应分析、综合分析
2.熟悉差分平稳序列的建模步骤。
3.掌握单位根检验、协整检验、动态回归模型的建立。
二、实验原理
1.序列的各种变化都归结于四大因素的综合影响：
长期趋势（Trend）,循环波动（Circle）,季节性变化（Season）,机波动（Immediate）.常假设它们有如下的相互模型：
加法模型
乘法模型
混合模型模型结构不唯一
2.非平稳序列如果能通过适当阶数的差分后实现平稳，就可以对差分后序列进行ARMA模型拟合了，所以ARIMA模型是差分运算与ARMA模型的组合
3.单位根检验：
（1）DF检验；
（2）ADF检验；（3）PP检验；
4.动态回归模型ARIMAX
如果两个非平稳序列之间具有协整关系，则先建立它们的回归模型，再对平稳的残差序列建立ARMA模型。
三、实验内容
1、P202页：
第7题（X11因素分解法）
2、P155页：
第3题（乘积季节模型）
3、P240页：
第4题
出口为
，进口为
，回答以下问题
（1）画出
，
的时序图，用单位根检验序列它们的平稳性；
（2）对
分别拟合模型（提示：
建立ARIMA模型）；
（3）考察
的协整关系，建立
的协整模型，同时建立误差修正模型。
四、实验过程
（一）P202页：
第7题（X11因素分解法）
1.绘制序列时序图。
（程序见附录）
由上图可得季节序列的振幅随序列水平的变化而变化，所以季节效应与趋势效应不独立，采用乘法模型:
2.进入x-11季节调整模型
经过三个阶段共十步的重复迭代后，得到如下的拟合效果图：
显然，该地区奶牛的月度产奶量序列具有显著的季节变动特征。
（二）P155页：
第3题（乘积季节模型）
1.绘制序列时序图。
绘制时序图，如图1所示（程序见附录1）。
图1美国月度事故死亡人数序列时序图
时序图显示该序列具有以年为周期的季节效应。
2.差分平稳化：
对原序列作1阶12步差分，希望提取原序列季节效应，差分后序列时序图如图2所示。
图2美国月度事故死亡人数1阶12步差分后序列时序图
时序图显示差分后序列类似平稳。
3.模型定阶：
考察差分后序列自相关图，如图3，进一步确定平稳性判断，并估计拟合模型的阶数。
图3美国月度事故死亡人数1阶12步差分后序列自相关图
自相关图显示延迟12阶自相关系数显著大于2倍标准差范围，这说明差分后序列中仍蕴含着非常显著的季节效应。
延迟1阶的自相关系数也大于2倍的标准差，这说明差分后序列还具有短期相关性。
观察偏自相关图，如图4，得到的结论和上面的结论一致。
图4美国月度事故死亡人数1阶12步差分后序列偏自相关图
图5序列白噪声检验
图5显示，原序列延迟各阶LB统计量的P值小于显著性水平0.05，所以拒绝原假设，序列不通过白噪声检验。
根据差分后序列的自相关图和偏自相关图的性质，拟合乘积季节模型
。
自相关图显示，12阶以内的自相关系数1阶截尾，偏自相关图显示，12阶以内的偏自相关系数1阶截尾，所以尝试使用ARMA（1，0）模型提取差分后序列的短期自相关信息。
再考虑季节自相关特征，这时考察延迟12阶、24阶等以周期长度为单位的自相关系数和偏自相关系数的特征。
自相关图显示延迟12阶自相关系数显著非零，而偏自相关图显示延迟12阶偏自相关系数显著非零，这时用以12步为周期的
模型提取差分后序列的季节自相关信息。
4.参数估计：
图6拟合模型
综合前面的差分信息，我们要拟合的乘积季节模型为
。
使用条件最小二乘估计方法，确定该模型的口径为：
5.模型检验：
对序列拟合
模型，模型及模型参数的显著性检验如图7、8所示。
图7模型参数的显著性
由图7知，拟合效果显示模型参数显著。
图8残差白噪声检验
对拟合模型进行白噪声检验，结果显示P值都大于显著性水平0.05.接受原假设，残差序列通过白噪声检验，模型显著，说明模型拟合良好，对序列相关信息提取充分。
将序列拟合值和序列观察值联合作图，如图9所示。
图9美国月度事故死亡人数拟合效果图
说明：
图中，点为序列观察值；曲线为序列拟合值。
从图9可以直观地看出该乘积季节模型对原序列的拟合效果良好。
（三）P240页：
第4题
1．画出
，
的时序图，用单位根检验序列的平稳性；
输出时序图如图1所示（程序见附录2）。
图1我国出口总额Xt、进口总额yt时序图
图1中，黑色为出口总额xt序列时序图，红色为进口总额yt序列时序图。
从图1中可以看出出口总额xt序列、进口总额yt序列均显著非平稳，这个直观判断还可以通过单位根检验验证。
同时时序图显示这两个序列具有某种同变关系。
对我国出口总额序列xt进行ADF检验，单位根检验结果如图2所示。
图2出口总额xt白噪声、单位根检验
检验结果显示，无论考虑何种类型的模型，检验统计量的P值均显著大于0.05的显著性水平，所以可以认为中国我国出口总额序列xt显著非平稳，且这六种处理均不能实现残差序列平稳。
对我国进口总额序列yt进行ADF检验，单位根检验结果如图3所示。
图3进口总额yt白噪声、单位根检验
同出口序列xt的检验结果一样，在显著性水平取为0.05时，可以认为我国进口序列yt非平稳，且这六种处理均不能实现残差序列平稳。
显然，这两个序列的ADF检验结果与根据时序图得到的直观判断完全一致
2．对
分别拟合模型（提示：
建立ARIMA模型）；
对我国出口对数序列lnxt和进口对数序列lnyt绘制时序图，如图4所示。
图4我国出口总额Xt、进口总额yt取对数时序图
图4中，黑色线代表我国出口对数序列lnxt，红色线代表我国进口对数序列lnyt。
时序图显示这两个对数序列有显著的上升趋势，为典型的非平稳序列。
同时时序图显示这两个序列具有某种同变关系。
因为序列呈现出近似线性趋势，所以选择1阶差分。
1阶差分后出口对数序列lnxt时序图如图5所示。
图5对数序列Lnx差分时序图
时序图显示，lnxt差分后序列在均值附近比较稳定地波动。
为了进一步确定平稳性，考察差分后序列的自相关图，如图6所示。
图6对数序列Lnxt差分后自相关图
自相关图显示序列有很强的短期相关性，所以可以初步认为lnxt1阶差分后序列平稳。
对平稳的1阶差分序列进行白噪声检验，白噪声检验结果如图7所示。
图7lnxt一阶差分后序列白噪声检验
在检验的显著性水平取为0.05的条件下，由于延迟6阶、12阶的P值均小于0.05，所以lnxt差分后的序列不能视为白噪声序列，即差分后序列还蕴含着不容忽视的相关信息可以提取。
对平稳非白噪声差分序列拟合ARMA模型，1阶差分后序列的自相关图（见图6）已经显示该序列有不截尾的性质。
再考察其偏自相关系数的性质，如图8所示。
图8对数序列Lnxt差分后偏自相关图
偏自相关图显示出1阶截尾性，所以考虑用AR
（1）模型拟合lnxt1阶差分后序列。
考虑到前面已经进行的1阶差分运算，实际上是用
模型拟合原序列。
对序列拟合
模型，模型参数及模型的显著性检验如图9、10所示。
图9模型参数显著性检验
由图9知，系数显著性检验显示两参数均显著。
对残差序列进行白噪声检验，检验结果如图10所示。
图10残差白噪声检验
显然，拟合检验统计量的P值都显著大于显著性检验水平0.05，可以认为残差序列即为白噪声序列，模型显著，这说明
模型对lnxt序列建模成功。
图11模型
在条件最小二乘估计原理下，拟合结果为：
2
将对数序列拟合值lnxt和对数序列观察值lnxt联合作图，如图12所示。
图12对数序列Lnxt拟合效果图
说明：
图中，星号为序列观察值；曲线为拟合值。
从图可以直观地看出该
模型对原序列的拟合效果良好。
因为对数序列lnyt呈现出近似线性趋势，所以选择1阶差分。
1阶差分后进口对数序列lnyt时序图如图13所示。
图13对数序列Lny差分时序图
时序图显示，lnyt差分后序列在均值附近比较稳定地波动。
为了进一步确定平稳性，考察差分后序列的自相关图，如图14所示。
图14对数序列Lnyt差分后自相关图
自相关图显示序列有很强的短期相关性，所以可以初步认为lnyt1阶差分后序列平稳。
对平稳的1阶差分序列进行白噪声检验，白噪声检验结果如图15所示。
图15lnyt一阶差分后序列白噪声检验
在检验的显著性水平取为0.05的条件下，由于延迟6阶的P值小于0.05，所以lnyt差分后的序列不能视为白噪声序列，即差分后序列还蕴含着不容忽视的相关信息可以提取。
对平稳非白噪声差分序列拟合ARMA模型，1阶差分后序列的自相关图（见图14）已经显示该序列有1阶截尾的性质。
再考察其偏自相关系数的性质，如图16所示。
图16对数序列Lnyt差分后偏自相关图
偏自相关图显示该序列1阶截尾的性质，所以考虑用AR
（1）模型拟合lnyt1阶差分后序列。
考虑到前面已经进行的1阶差分运算，实际上是用
模型拟合原序列。
对序列拟合
模型，模型参数及模型的显著性检验如图17、18所示。
图17模型参数显著性检验
由图17知，系数显著性检验显示两参数均显著。
对残差序列进行白噪声检验，检验结果如图18所示。
图18残差白噪声检验
显然，拟合检验统计量的P值都显著大于显著性检验水平0.05，可以认为残差序列即为白噪声序列，模型显著。
这说明
模型对该序列建模成功。
图19模型
在条件最小二乘估计原理下，拟合结果为：
将对数序列拟合值lnyt和对数序列观察值lnyt联合作图，如图20所示。
图20对数序列Lnyt拟合效果图
说明：
图中，星号为序列观察值；曲线为拟合值。
从图20可以直观地看出该
模型对原序列的拟合效果良好。
3．考察
的协整关系，建立
的协整模型，同时建立误差修正模型。
对我国出口对数序列lnxt和进口对数序lnyt绘制时序图，如图4所示。
可以发现时序图显示这两个序列具有某种同变关系，可以考虑建立ARIMAX模型。
对lnxt、lnyt、lnxt1阶差分、
序列分别进行单位根检验（ADF）。
输出结果如图21——24所示。
图21对数序列lnXt1阶单位根检验
检验结果显示，无论考虑何种类型的模型，检验统计量的P值均显著大于0.05的显著性水平，所以可以认为中国我国出口总额对数序列lnxt显著非平稳，且这六种处理均不能实现残差序列平稳。
图22对数序列lnyt1阶单位根检验
同出口对数序列lnxt的检验结果一样，在显著性水平取为0.05时，可以认为我国进口对数序列lnyt非平稳，且这六种处理均不能实现残差序列平稳。
图23
1阶单位根检验
检验结果显示，无论考虑何种类型的模型，检验统计量的P值均显著小于0.05的显著性水平，拒绝原假设，所以可以认为中国我国出口总额
对序列显著平稳。
图24
1阶单位根检验
同出口
序列的检验结果一样，在显著性水平取为0.05时，可以认为我国进口
序列平稳，且这六种处理均能实现残差序列平稳。
利用最小二乘估计，回归模型输出结果如图25所示。
图25回归模型结果
构造出的回归模型如下：
图26对数序列Lnyt与对数序列lnxt之间的相关图
相关图显示对数序列lnyt在延迟阶数为零时与对数序列lnxt相关关系最大。
因此可以将对数序列lnyt与对数序列lnxt同期建模。
图27残差单位根检验
残差序列平稳，说明对数序列lnyt与对数序列lnxt之间具有协整的关系，我们可以大胆的在这两个对数序列之间建立动态回归模型而不必担心虚假回归问题。
图28模型参数显著性检验，无常数项
考察残差序列白噪声检验结果，如图29所示。
图29残差序列白噪声检验
输出结果显示，延迟各阶LB统计量的P值都大于显著性水平0.05，可以认为残差序列为白噪声序列，对序列相关信息提取充分。
图30模型
根据输出的模型拟合结果可知，最后的拟合模型口径为：
将对数序列拟合值lnyt和对数序列观察值lnyt联合作图，如图31所示。
图31对数序列Lnyt拟合效果图
说明：
图中，星号为序列观察值；曲线为拟合值。
从图可以直观地看出该ARIMAX模型对原序列的拟合效果良好。
将序列拟合值yt和序列观察值yt联合作图，如图32所示。
图32yt拟合效果图
说明：
图中，星号为序列观察值；曲线为拟合值。
从图可以直观地看出该ARIMAX模型对原序列的拟合效果良好。
构造误差修正模型:
为了研究我国的进出口总额的短期波动特征，我们利用差分序列
和
以及前期误差序列
，构造ECM模型：
输出结果如图33所。
图33误差修正模型显著性检验
由图33输出结果结果知ECM模型为：
附录
程序1：
dataexample4_7;
inputx@@;
t=intnx（'quarter','1jan1978'd,_n_-1）;
formattyyq4.;
cards;
589561640656727697640599568577553582
600566653673742716660617583587565598
628618688705770736678639604611594634
658622709722782756702653615621602635
677635736755811798735697661667645688
713667762784837817767722681687660698
717696775796858826783740701706677711
734690785805871845801764725723690734
750707807824886859819783740747711751
;
procx11data=example4_7;
quarterlydate=t;
varx;
outputout=outb1=xd10=seasond11=adjustedd12=trendd13=irr;
dataout;
setout;
estimate=trend*season/100;
procgplotdata=out;
plotseason*t=2adjusted*t=2trend*t=2irr*t=2;
plotx*t=1estimate*t=2/overlay;
symbol1c=blacki=joinv=star;
symbol2c=redi=joinv=nonew=2;
run;
程序2：
datati5_3;
inputx@@;
dif1_12=dif12（dif（x））;
time=intnx（'month','1jan1973'd,_n_-1）;
formattimeyear4.;
cards;
9007.008106.008928.009137.0010017.0010826.0011317.0010744.00
9713.009938.009161.008927.007750.006981.008038.008422.008714.009512.00
10120.009823.008743.009129.008710.008680.008162.007306.008124.00
7870.009387.009556.0010093.009620.008285.008433.008160.008034.00
7717.007461.007776.007925.008634.008945.0010078.009179.008037.00
8488.007874.008647.007792.006957.007726.008106.008890.009299.0010625.00
9302.008314.008850.008265.008796.007836.006892.007791.008129.009115.00
9434.0010484.009827.009110.009070.008633.009240.00
;run;
procgplot;
plotx*time=1dif1_12*time=2;
symbol1c=coralv=circlei=join;
symbol2c=bluev=stari=join;run;
procarima;
identifyvar=x（1,12）;
estimatep=1q=
（1）（12）;
forecastlead=0id=timeout=out;run;
procgplotdata=out;
plotx*time=1forecast*time=2/overlay;
symbol1c=blacki=nonev=doth=0.2;
symbol2c=redi=joinv=none;
run;
程序3:
datati6_4;
inputyearxtyt;
lnxt=log（xt）;
lnyt=log（yt）;
diflnx=dif（lnxt）;
diflny=dif（lnyt）;
cards;
19502021.3
195124.235.3
195227.137.5
195334.846.1
19544044.7
195548.761.1
195655.753
195754.550
19586761.7
195978.171.2
196063.365.1
196147.743
196247.133.8
19635035.7
196455.442.1
196563.155.3
19666661.1
196758.853.4
196857.650.9
196959.847.2
197056.856.1
197168.552.4
197282.964
1973116.9103.6
1974139.4152.8
1975143147.4
1976134.8129.3
1977139.7132.8
1978167.6187.4
1979211.7242.9
1980271.2298.8
1981367.6367.7
1982413.8357.5
1983438.3421.8
1984580.5620.5
1985808.91257.8
19861082.11498.3
198714701614.2
19881766.72055.1
198919562199.9
19902985.82574.3
19913827.13398.7
19924676.34443.3
19935284.85986.2
199410421.89960.1
199512451.811048.1
199612576.411557.4
199715160.711806.5
199815223.611626.1
199916159.813736.5
200020634.418638.8
200122024.420159.2
200226947.924430.3
200336287.934195.6
200449103.346435.8
200562648.154273.7
200677594.663376.9
200793455.673284.6
2008100394.979526.5
;run;
procgplotdata=ti6_4;
plotxt*year=1yt*year=2/overlay;
symbol1c=blacki=joinv=none;
symbol2c=redi=joinv=nonew=2l=2;run;
procarimadata=ti6_4;
identifyvar=xtstationarity=（adf=1）;
identifyvar=ytstationarity=（adf=1）;run;
procgplotdata=ti6_4;
plotlnxt*year=1lnyt*year=2/overlay;
plotdiflnx*year=1diflny*year=2;
symbol1c=blacki=joinv=circle;
symbol2c=redi=joinv=star;run;
procarima;
identifyvar=lnxt
（1）;
estimatep=1;
forecastlead=0id=yearout=out1;
identifyvar=lnyt
（1）;
estimatep=1;
forecastlead=0id=yearout=out2;run;
procgplotdata=out1;
plotlnxt*year=1forecast*year=2/overlay;
symbol1c=blacki=nonev=star;
symbol2c=redi=joinv=none;run;
procgplotdata=out2;
plotlnyt*year=1forecast*year=2/overlay;
symbol1c=blacki=nonev=star;
symbol2c=redi=joinv=none;run;
procarimadata=ti6_4;
identifyvar=lnxtstationarity=（adf=1）;
identif
升级会员 