時間序列整合分析(一)

1、1. 缺省值的補足：2. 時序圖：（檢驗平穩(wěn)性）3. 自相關(guān)函數(shù)：（檢驗平穩(wěn)性）4. 計算標準正態(tài)分布的概率：5. 計算標準正態(tài)分布的分位數(shù)：6. 計算標準t分布的概率7. 計算標準t分布的分位數(shù)8. 計算標準F分布的概率9. 計算標準F分布的分位數(shù)10. 計算標準卡方分布的概率11. 計算標準卡方分布的分位數(shù)12. 方差的同齊性檢驗：將數(shù)據(jù)進行適當分組，這里將4個分為一組，一共四組PrF的值大于0.05 故接受H0，認為各組方差之間沒有顯著的差異。13. 方差的同質(zhì)性檢驗：將數(shù)據(jù)進行適當分組，這里將4個分為一組，一共四組 根據(jù)上結(jié)果列出方差分析表： 方差來源平方和自由度均方和F值顯著性A誤差

2、878376140372932029279270186.44.17總和228210523 F的p值小于0.05 我們認為原始數(shù)據(jù)方差不同質(zhì)。14. 序列的白噪聲檢驗（檢驗純隨機性）： 可以看出，LB(6)=95.84，其p值小于0.05； LB(12)=190.40，其p值小于0.05；顯然該序列不是白噪聲序列，即不是純隨機性序列。（p值都大于0.05時才是純隨機序列）15. 平穩(wěn)序列的自相關(guān)函數(shù)和偏自相關(guān)函數(shù)的形式：（沒有程序的） 模型AR(p)MA(q)ARMA(p,q)ACF自相關(guān)拖尾截尾拖尾PACF偏自相關(guān)截尾拖尾拖尾16. 一個例子：（利用平穩(wěn)序列建模 進行預(yù)測） 我國1975-20

3、06年GDP的年增長率為下表（數(shù)據(jù)略），對我國1975-2006年GDP的年增長率進行建模，并對2007至2011年我國的GDP增長率進行預(yù)測。 (1)首先畫出我國1975-2006年GDP增長率的時序圖。data ex;input x;t=_n_;cards;8.7 -1.6 7.6 11.7 7.6 7.8 5.2 9.1 10.9 15.2 13.5 8.8 11.6 11.3 4.1 3.8 9.2 14.2 14 13.1 10.9 10 9.3 7.8 7.6 8.4 8.3 9.1 10 10.1 10.4 10.7;proc gplot;symbol i=jiont v=dot

4、;plot x*t;run;從圖中直觀的可以看出有奇異點（2） 將奇異點看成缺省值，利用以下程序來求缺省點的值：data ex;input x;time=intnx(month,01jan1975d,_n_-1);format time data;cards;8.7 . 7.6 11.7 7.6 7.8 5.2 9.1 10.9 15.2 13.5 8.8 11.6 11.3 4.1 3.8 9.2 14.2 14 13.1 10.9 10 9.3 7.8 7.6 8.4 8.3 9.1 10 10.1 10.4 10.7;proc expand data=ex out=ex1;id time

5、;proc print data=ex1;run;結(jié)果可知，缺省點的值為2.4（3） 利用修正后的數(shù)據(jù)再進行時序分析，根據(jù)以下程序：可以看出GDP增長率修正后的數(shù)據(jù)序列平穩(wěn)。BIC（5,0）=-0.24488的值最小，考慮建立AR(5)模型。（4） 模型的建立data ex;input x;time=intnx(month,01jan1975d,_n_-1);format time year4.;cards;8.7 2.4 7.6 11.7 7.6 7.8 5.2 9.1 10.9 15.2 13.5 8.8 11.6 11.3 4.1 3.8 9.2 14.2 14 13.1 10.9 10

6、 9.3 7.8 7.6 8.4 8.3 9.1 10 10.1 10.4 10.7;proc arima;identify var=x nlag=12 minic p=(0:5) q=(0:5);estimate p=5;run;從上圖中可以看出，有些參數(shù)不顯著，我們將其去掉，建立最精干的模型。（其中可以看出，AR1,3 AR1,4 AR1,5 的p值遠遠大于0.05）所以，將estimate p=5改為estimate p=(1,2)，即程序為：data ex;input x;time=intnx(month,01jan1975d,_n_-1);format time year4.;car

7、ds;8.7 2.4 7.6 11.7 7.6 7.8 5.2 9.1 10.9 15.2 13.5 8.8 11.6 11.3 4.1 3.8 9.2 14.2 14 13.1 10.9 10 9.3 7.8 7.6 8.4 8.3 9.1 10 10.1 10.4 10.7;proc arima;identify var=x nlag=12 minic p=(0:5) q=(0:5);estimate p=(1,2);run;可見所有的p 值都小于0.05 通過了檢驗。所以模型為：最后AR（5）模型的殘差分析為（即模型的殘差白噪聲檢驗）：LB(6)=3.45 其p值為0.4851大于0.0

8、5 ，故通過檢驗（其他的也是類似），所以該模型的擬合效果很好。（5） 用此模型做預(yù)測data ex;input x;t=intnx(year,01jan1975d,_n_-1);format t year4.;cards;8.7 2.4 7.6 11.7 7.6 7.8 5.2 9.1 10.9 15.2 13.5 8.8 11.6 11.3 4.1 3.8 9.2 14.2 14 13.1 10.9 10 9.3 7.8 7.6 8.4 8.3 9.1 10 10.1 10.4 10.7;proc arima;identify var=x nlag=12 minic p=(0 5) q=(0

9、 5);estimate p=(1,2) method=cls;forecast lead=5 id=t out=results;proc gplot data=results;plot x*t=1 forecast*t=2 l95*t=3/overlay;symbol1 c=blue i=jion v=star;symbol2 c=red i=jion v=none l=1 w=1;symbol3 c=green i=jion v=none l=2 w=2;run;可以看出，原始數(shù)據(jù)絕大部分都在預(yù)測區(qū)域內(nèi)，而且越是近期的數(shù)據(jù)離預(yù)測曲線越近，表明模型建立的比較合理，預(yù)測效果比較精準。17. 線

10、性擬合：序列呈現(xiàn)出顯著的線性特征我國1995-2006年期刊種類數(shù)據(jù)如表（數(shù)據(jù)略），試建模并進行預(yù)測。（1） 畫出改時間序列的散點圖data ex;input x;t=_n_;cards;7583 7916 7918 7999 8187 8725 8889 9029 9074 9490 9468 9468;proc gplot;symbol i=jiont v=dot;plot x*t;run;發(fā)現(xiàn)該序列有顯著的線性遞增趨勢，于是考慮使用線性模型：來擬合該趨勢的發(fā)展。（2） 擬合模型，進行參數(shù)估計和檢驗。data ex;input x;t=_n_;cards;7583 7916 7918 79

11、99 8187 8725 8889 9029 9074 9490 9468 9468;proc reg;model x=t;run;在t檢驗中，p值都小于0.05 拒絕H0 接受H1，認為參數(shù)都顯著非0。在F檢驗中，p值小于0.05 認為方程顯著。得到趨勢模型：（3） 對殘差進行白噪聲檢驗。data ex;input x;t=_n_;epsilon=x-7423.40909-188.01399*t;cards;7583 7916 7918 7999 8187 8725 8889 9029 9074 9490 9468 9468;proc arima;identify var=epsilon n

12、lag=12;run;LB(6)=5.88 其p值為0.4367，遠遠大于0.05，故通過參差檢驗，殘差為白噪聲序列。表明模型信息提取充分。我們將利用前面得到的趨勢模型進行預(yù)測，將t=13 t=14分別代入模型中，可得：2007年我國期刊種類預(yù)測值：7423.40909+188.01399*13=9867.592008年我國期刊種類預(yù)測值：7423.40909+188.01399*14=10055.6018. 曲線擬合：序列在發(fā)展過程中呈現(xiàn)出以不同的速率增長或下降，或者由逐漸增長到逐漸衰退等各種不同的性態(tài)，即非線性特征對1980-2007年GDP的數(shù)據(jù)進行建模和預(yù)測（數(shù)據(jù)略）（1）畫出1980

13、-2007年全國GDP散點圖：data ex;input GDP;t=intnx(year,01jan1980d,_n_-1);format t year4.;cards;4517.8 4862.4 5294.7 5934.5 7171 8964.4 10202.211962.5 14928.3 16909.2 18547.9 21617.8 26638.1 34634.446759.4 58478.1 67884.6 74462.6 78345.2 82067.5 89442.297314.8 118020.7 135822.8 159878.3 183217.4 211923 249530

14、;proc gplot;plot GDP*t;symbol c=blue i=jion v=dot;run;由時序圖可以看出序列呈現(xiàn)指數(shù)增長趨勢，因此，我們提出初步模型（2） 用最小二乘法進行參數(shù)估計：data ex;input x_t;t=_n_;cards;4517.8 4862.4 5294.7 5934.5 7171 8964.4 10202.211962.5 14928.3 16909.2 18547.9 21617.8 26638.1 34634.446759.4 58478.1 67884.6 74462.6 78345.2 82067.5 89442.297314.8 1180

15、20.7 135822.8 159878.3 183217.4 211923 249530;proc nlin method=gauss;model x_t=a+b*exp(c*t);parameters a=1 b=0.5 c=0.5;der.a=1;der.b=exp(c*t);der.c=t*b*exp(c*t);output predicted=x_that out=result;proc gplot data=result;plot x_t*t=1 x_that*t=2/overlay;symbol1 c=black i=none v=dot l=1 w=1;symbol2 c=bl

16、ue i=jion v=dot l=2 w=2;Run;得到非線性模型的擬合圖：非線性模型的參數(shù)估計通過了檢驗。（3） 對殘差進行白噪聲檢驗：data ex;input x_t;t=_n_;epsilon=x_t+3989.1-7252.7*exp(0.1305*t);cards;4517.8 4862.4 5294.7 5934.5 7171 8964.4 10202.211962.5 14928.3 16909.2 18547.9 21617.8 26638.1 34634.446759.4 58478.1 67884.6 74462.6 78345.2 82067.5 89442.297

17、314.8 118020.7 135822.8 159878.3 183217.4 211923 249530;proc arima;identify var=epsilon nlag=12;run;可以看出，LB(6)=54.08其p值小于0.05，沒有通過檢驗。 LB(12)=60.12其p值小于0.05，沒有通過檢驗。表明模型參差不是白噪聲序列，信息提取不充分，我們還需要對模型進行進一步改進。（4） 由（1）中的散點圖可以看出，序列具有一定的指數(shù)趨勢，首先對其取對數(shù)：data ex;input GDP;y=log(GDP);t=intnx(year,01jan1980d,_n_-1);f

18、ormat t year4.;cards;4517.8 4862.4 5294.7 5934.5 7171 8964.4 10202.211962.5 14928.3 16909.2 18547.9 21617.8 26638.1 34634.446759.4 58478.1 67884.6 74462.6 78345.2 82067.5 89442.297314.8 118020.7 135822.8 159878.3 183217.4 211923 249530;proc gplot;plot y*t;symbol c=blue i=jion v=dot;run;取對數(shù)后的序列呈現(xiàn)線性增長

19、趨勢，對取對數(shù)后的數(shù)據(jù)進行一階差分：data ex;input GDP;y=log(GDP);y1=dif(y);t=intnx(year,01jan1980d,_n_-1);format t year4.;cards;4517.8 4862.4 5294.7 5934.5 7171 8964.4 10202.211962.5 14928.3 16909.2 18547.9 21617.8 26638.1 34634.446759.4 58478.1 67884.6 74462.6 78345.2 82067.5 89442.297314.8 118020.7 135822.8 159878.

20、3 183217.4 211923 249530;proc gplot;plot y1*t;symbol c=blue i=jion v=dot;run;取對數(shù)后為平穩(wěn)序列，接下來按照16中平穩(wěn)序列的處理方法進行建模和預(yù)測。19. 季節(jié)指數(shù)水平法（適用于無明顯上升或下降變動趨勢，主要受季節(jié)變動和不規(guī)則變動影響的 時間序列）下面是臨海大藥房季節(jié)藥博利康尼2004-2007年各季度的銷售量。 季度2004200520062007190110103982605666583253532364120155123115（1） 先做出時序圖：data ex;input x;t=_n_;format t ye

21、ar4.;cards;90 60 25 120 110 56 35 155 103 66 32 123 98 58 36 115;proc gplot;plot x*t;symbol c=blue i=jion v=dot;run;由上圖可以看出，臨海大藥房季節(jié)藥博利康尼為季節(jié)性藥品，在不同的季度其銷售量的波動很明顯，出現(xiàn)明顯的銷售旺季和銷售淡季，改時間序列為非平穩(wěn)的。（2） 得到藥品銷售各平均值和季節(jié)指數(shù)如下： 總平均=80.125 第一季度平均=100.25 第二季度平均=60.00 第三季度平均=32.00 第四季度平均=128.25季節(jié)指數(shù)S1=100.25/80.125=1.2511

22、7 S2=0.748830 S3=0.399376 S4=1.600624可以看出季節(jié)指數(shù)差別很大，說明藥品銷售量有明顯的季節(jié)性。（3） 對2008年進行預(yù)測：進行預(yù)測事，一般選擇最近年份的平均值，這里，2007年的每季度的平均值為（98+58+36+115）/4=76.75第一季度：76.75*1.25117=95.9375第二季度：76.75*0.74883=57.5625第三季度：76.75*0.399376=30.7第四季度：76.75*1.600624=122.8全年：95.9375+57.5625+30.7+122.8=30720. 季節(jié)指數(shù)趨勢法：如果序列既存在季節(jié)變動，同時各年

23、水平或同月（同季）水平呈現(xiàn)上升或下降趨勢，這是應(yīng)該采用季節(jié)指數(shù)趨勢法。K期的移動平均：K期中心移動平均： 最后計算：例：下表為某市某品牌桶裝水銷售記錄數(shù)據(jù)，我們據(jù)此數(shù)據(jù)做出2006年該桶裝水的銷售記錄的預(yù)測。年份季度銷售量Xt/萬桶4期中心化移動平均季節(jié)比率消除季節(jié)影響的銷售量Xt/St200114659.4595126358.1637738862.51.40864.6961245163.8750.7984364.38623200215065.1250.7677564.629927066.251.056664.6264139168.251.333366.9016745770.50.808517

24、1.96108200316072.50.8275977.5558827874.251.0505172.0122939974.6251.3266372.7831346375.6250.8330679.53593200415778.3750.7272773.6780928979.3751.1212682.16787311080.1251.3728580.8701546082.3750.7283875.74851200516685.50.7719385.3114729888.751.1042390.47698312692.6330847088.37326（1） 先做出序列的散點圖：data ex;i

25、nput x;t=intnx(quarter,01jan2001d,_n_-1);format t yyq4.;cards;466388515070915760789963578911060669812670;proc gplot;plot x*t;symbol c=blue i=jion v=dot;run;觀察該市桶裝水銷售量序列由長期趨勢和一年為周期的季節(jié)波動共同影響，因此要考慮剔除趨勢再分析季節(jié)特征。（2） 剔除趨勢：NO1：以一年的季度數(shù)4為k，對時間序列觀測值進行k項中心移動平均，形成新的序列，如表中第四列所示，以此為長期趨勢Tt。在消除了時間序列的季節(jié)和不規(guī)則性之后，可以看出，具

26、有很明顯的趨勢性。 data ex;input x; t=intnx(quarter,01jan2001d,_n_-1);format t yyq4.; cards; 62.563.87565.12566.2568.2570.572.574.2574.625 75.62578.37579.37580.12582.37585.588.75 ; proc gplot;plot x*t;symbol c=blue i=jion v=dot;run;NO2：將各期觀察值除以同期移動均值作為季節(jié)比率，即表中倒數(shù)第二列。NO3：各年同季的季節(jié)比率平均，季節(jié)平均比率可以消除不規(guī)則變動，即得其季節(jié)指數(shù)St。

27、做法與19中（2）中的做法是一樣的。 得到：季度季節(jié)不規(guī)則成分的數(shù)值（StIt）季節(jié)指數(shù)（St）10.76775 0.82759 0.72727 0.771930.7736421.0566 1.05051 1.12126 1.104231.0831531.408 1.3333 1.32663 1.372851.3602140.79843 0.80851 0.83306 0.728380.79209季節(jié)指數(shù)模型需要平均季節(jié)指數(shù)等于1.00，在表中，季節(jié)指數(shù)總和為4.009084，平均值基本等于1，所以不用調(diào)整。對于不等于1的，要進行一些調(diào)整：用每一個季節(jié)指數(shù)乘以季度總和再除調(diào)整之前的季節(jié)指數(shù)之和

28、。NO4：用時間序列的每個觀察值處以相應(yīng)的季節(jié)指數(shù)，消除時間序列季節(jié)影響，得到趨勢序列，即，表中最后一列。（3） 利用消除季節(jié)影響后的線性趨勢預(yù)測值建立回歸模型，進行回歸分析與白噪聲檢驗。data ex;input T_t;t=_n_;format t yyq4.;cards;59.4595158.1637764.6961264.3862364.629964.6264166.9016771.9610877.5558872.0122972.7831379.5359373.6780982.1678780.8701575.7485185.3114790.4769892.6330888.37326;p

29、roc reg;model T_t=t;run;趨勢模型的參數(shù)和模型顯著，都通過了檢驗。其表達式為：進行白噪聲檢驗：data ex;input T_t;t=_n_;format t yyq4.;res=T_t-57.05063-1.64266*t;cards;59.4595158.1637764.6961264.3862364.629964.6264166.9016771.9610877.5558872.0122972.7831379.5359373.6780982.1678780.8701575.7485185.3114790.4769892.6330888.37326;proc arima

30、;identify var=res nlag=12;run;可以看出，白噪聲檢驗通過。趨勢模型是適應(yīng)的。（4） 進行預(yù)測和季節(jié)調(diào)整。對同時有趨勢和季節(jié)成分的時間序列，進行預(yù)測的最后是用季節(jié)指數(shù)調(diào)整趨勢預(yù)測值。2006.1：Tt=57.05063+1.64266*17=91.54649 調(diào)整后：91.54649*0.77364=70.82362006.2：Tt=57.05063+1.64266*18=93.18915 調(diào)整后：93.18915*1.08315=100.9同理2006.3 2006.4。21. X11方法（變量的第二類調(diào)查統(tǒng)計方法季節(jié)調(diào)整方案）：比季節(jié)指數(shù)法預(yù)測更精確例：某市醫(yī)院三

31、年中各個季度接受腫瘤治療的人次如表所示，是進行季節(jié)調(diào)整，計算出季節(jié)指數(shù)，并進行預(yù)報。季度1234199920002001142451461014976143351470115066144261479215066145181488415249（1）data ex;input x;t=intnx(quarter,01jan1999d,_n_-1);format t yyq.4;cards;142451433514426145181461014701147921488414976150661506615249;proc x11 data= ex; quarterly date=t;output out=out b1=x d10=season d11=adjusted d12=trend d13=irr;/*b1為原序列d10列出季節(jié)指數(shù)d11列出季節(jié)調(diào)整后的序列 d12趨勢擬合 d13列出最后的不規(guī)則波動值*/run;圖1 生成的原始數(shù)據(jù)表圖2 時間序列的季節(jié)因子St（%）圖3 季節(jié)調(diào)整后的序