时间序列回归模型1干预分析1.1概念及模型Box和Tiao引入的干预分析提供了对于干预影响时间序列的效果进行评估的一个框架,假设干预是可以通过时间序列的均值函数或者趋势而对过程施加影响,干预可以自然产生也可以人为施加的,如国家的宏观调控等。
其模型可以如下表示:其中mt代表均值的变化,Nt是ARIMA过程。
1.2干预的分类阶梯响应干预區案1“ 書聲新镖第应干严的苕爭第见複也[榔帝右一牛时闽单恆的延遇)01 "4》* a_e—4 f-辜—右4—*—T1)诅畠严to it r ■P■1FV*1脉冲响应干预图聲1J4荷关脉冲愉血于预的一牲常见棋型(都带衬一个时伺单也的延迟)1.3干预的实例分析1.3.1 模型初探对数化航空客运里程的干预模型的估计现任回到每月航空客运蚩程的数据.如前所述’ 2(X)1年9刀的悉怖裳击事杵便航空客运徘徊于萧条之中,该T•预效应可用在200]年9月有脉亦输入的AR (1)过程柬表示*这一意外爭件对航克容运虽即时造底了一种强烈的激冷效应*因此*对此干预效应<9-11 »应)建模如下’叭=咖戶汙十1 3'严1 —M M展中,T代表2001年9小在这一衷示中*纽+助代表即时的9/11效应・且当^>1时* 纳(毗尸代表9门1效应对苴后A个月粉所造成的影响.这里还需要确定華础无扰过思的季节ARTMA 构*基于预干预数据,輛用一个AR1MA (0, 1, l)X<0・1, 0儿模型表示未愛扰的过程I券见图表11-5<> data(airmiles)>acf(as.vector(diff(diff(wi ndow(log(airmiles),e nd=c(2001,8)),12))),lag.max=48)# 用window 得到在911事件以前的未爱干预的时间序列子集Seri»es碍皿伽〔aimaiffi(響¥蹄[嚅律「皿"河,enc, =口起M 刖人对暂用的模型进行诊断>fitmode<-arima(airmiles,order=c(0,1,1),seas onal=list(order=c(0,1,0)))> tsdiag(fitmode)1.3.2拟合带有干预信息的模型函数:arimax(x, order = c(0, 0, 0), seas onal = list(order = c(0, 0, 0), period =NA),xreg = NULL, i nclude.mea n = TRUE, tran sform.pars = TRUE, fixed = NULL, in it = NULL, method = c("CSS-ML", "ML", "CSS"), n.cond, optim.c ontrol = list(), kappa = 1e+06, io = NULL, xtra nsf, tran sfer = NULL)arimax 函数扩展了 arima 函数,可以处理时间序列中干扰分析及异常值。
假设干扰影响 过程的均值,相对未受干扰的无价值函数的偏离用一些协变量的ARMA 滤波器的输出这种来 表示,偏差被称作传递函数。
构造传递函数的协变量通过xtransf参数以矩阵或者data.frame 的形式代入 arimax 函数。
air.m1=arimax(log(airmiles),order=c(0,1,1),seas on al=list(order=c( 0,1,1),period=12),xtra nsf=data.frame(l911=1*(seq(airmiles)==69), I911=1*(seq(airmiles)==69)),tran sfer=list(c(0,0),c(1,0)),xreg=data.frame(Dec96=1*(seq(airmile s)==12),Jan 97=1*(seq (ai rmiles)==13),Dec02=1*(seq(airmiles)==84)),method=' ML') > air.m1 Call:从诊断图可以看出存在三个异常点, acf 在12阶存在高度相关因此在季节中加入 MA (1)系数。
Coefficie nts: ma1 sma1 Dec96 Jan 97 Dec02 I911-MA0 I911.1-AR1 I911.1-MA0-0.3825 -0.6499 0.0989 -0.0690 0.0810-0.09490.8139-0.2715s.e. 0.0926 0.1189 0.0228 0.0218 0.02020.04620.09780.0439sigma A 2 estimated as 0.0006721: log-423.98画图plot(log(airmiles),ylab="log(airmiles)") poi nts(fitted(air.m1))Nin e11p=1*(seq(airmiles)==69) plot(ts(Ni ne11p*(-0.0949)+filter(Ni ne11p,filter=.8139,method='recursive',side=1)*(-0.2715), freque ncy=12,start=1996),type='h',ylab='9/11 Effects') abli ne(h=0)arimax(x = log(airmiles), order = c(0, 1, 1), seas onal = list(order =c(0, 1,1), period =12), xreg = data.frame(Dec96 = 1 * (seq(airmiles)12), Ja n97 = 1 * (seq(airmiles) == 13), Dec02 =1 * (seq(airmiles)==84)),="ML",xtra nsf = data.frame(I911 = 1 * (seq(airmiles)= :=69), I911 = 1* (seq(airmiles)==69)), tran sfer = =list(c(0, 0), c(1,0)))methodlikelihood = 219.99, aic从上图可以看出在 2003年底后,911事件的影响效应才平息,航班客运量恢复了正常。
2异常值在时间序列中异常有两种,可加异常和新息异常,分别记AO 和10。
2.1异常值示例 2.1.1模拟数据模拟一般的ARIMA ( 1,0,1 ),然后故意将第10个观测值变成异常值 10.> set.seed(12345)Start = 1End = 100Freque ncy = 1[1] 0.49180881 -0.22323665 -0.99151270 -0.733878180.51869129 1.86210605 2.19935472 2.60210165[17] 0.79130003 0.26265426 2.93414857 3.99045889 3.60822678 1.17845765 -0.87682948 -1.20637799[25] -1.39501221 -0.18832171 1.22999827 1.46814850 2.66647491 3.23417469 2.60349624 1.49513215[33] 1.48852142 0.95739219 1.30011654 1.73444053 2.84825103 3.73214655 4.23579456 3.37049790[41] 2.02783955 1.41218929 -0.29974176 -1.58712591 -1.34080878 0.10747609 1.44651081 1.67809487啦1W 20002HH 200*Times L.o-n.o-rt9二济[57] 1.70668201 1.37518194 1.91824534 0.14254056-2.88169481 -3.30372327 -1.74068408 -3.24868057[73] 2.00559443 0.86443324 0.46847572 0.723384981.60215098 1.25922277 1.53180859 0.96289779[81] 1.07712188 1.42386354 0.56318008 -0.46689543 -0.91861106 -1.92947085 -2.18188785 -1.02759087[89] 2.31088272 3.13847319 3.01237881 3.434548072.31539494 2.44909873 2.91589141 1.12648908[97] -0.08123871 0.44412579 0.26116418 -0.45815484 > y[10]<-102.1.2 模型初步判断> acf(y)Series yLag> pacf(y)Strlfl* 丫> eacf(y) AR/MA1 o o o o o o o o o o o o o o2 o o o o o o o o o o oo o o3 o x o o o o o o o o oo o o4 o x o o o o o o o o oo o o5 x x o o o o o o o o oo o o6 x o o o o o o o o o oo o o7 o x o o o o o o o o oo o o从三个的结果来看,可以初步分析y是AR (1)模型2.1.3 对模型时行拟合> m1=arima(y,order=c(1,0,0))> m1Call:arima(x = y, order = c(1, 0, 0))Coefficie nts:ar1 in tercept0.5419 0.7096s.e. 0.0831 0.36032.1.4 对模拟模型进行异常值探测> detectAO(m1)[,1] [,2] [,3]ind 9.000000 10.000000 11.000000Iambda2 -4.018412 9.068982 -4.247367> detectAO(m1,robust=F)[,1] Iind 10.000000Iambda2 7.321709> detectlO(m1)[,1] [,2]ind 10.000000 11.00000lambda1 7.782013 -4.67421AO探测结果认为第9 , 10 , 11.可能出现异常值。