Embed Size (px)
Citation preview
Modelos de Séries Temporais no Stata∗
Henrique Dantas Neder - Professor Associado Universidade Federal de Uberlândia
1 Introdução a Modelos Arima
O comando arima do Stata estima modelos ARIMA padrões que são autoregres-sivos na variável dependente e modelos estruturais com perturbações (erros)ARMA. Considere um processo autoregressivo de médias móveis de primeiraordem. O comando arima estima todos os parâmetros no modelo:
yt = xtβ + µt equação estruturalµt = ρµt−1 + θεt−1 + εt perturbação aleatória (erro), ARMA(1,1)onde:ρé o parâmetro de autocorrelação de primeira ordemθé o parâmetro de médias móveis de primeira ordemεt ~i.i.d. N(0,σ2), signi�cando que εt é um ruido brancoPodemos combinar as duas equações e escrever um modelo na forma geral
ARMA(p,q) no processo de perturbações como:
yt= xt β+ρ1(yt−1 − xt−1β)+ρ2(yt−2 − xt−2β) + ...+ ρp(yt−p − xt−pβ) (1)
+θ1εt−1 + θ2εt−2 + ...+ θqεt−q + εt
É também comum escrever a forma geral do modelo ARMA de uma formasucinta usando a notação do operador de lag:
ρ(Lp)(yt − xtβ) = θ(Lq)εtonde:Lp = 1− ρ1L− ρ2L2 − ...− ρpLpLq = 1 + θ1L+ θ2L
2 + ...+ θqLp
e Ljyt = yt−j
2 Modelos ARIMA
Modelos ARIMA puros sem a componente estrutural não tem regressores e sãofrequentemente escritos como autoregressões na variável dependente ao invés de
∗Este texto é baseado em uma tradução livre de diversas partes do Manual de SériesTemporais do Stata (que acompanha o software e pode ser visualizado através do menu Help- Documentation)
1
autoregressões na perturbação (erro) aleatório de uma equação estrutural. Porexemplo, um modelo ARMA(1,1) pode ser escrito como:
yt = α+ ρyt−1 + θεt−1 + εtExemplo 1: Modelo ARIMAEnders (2004, 87�93) considera um modelo ARMA para a série de preços
ao atacado dos Estados Unidos (WPI) usando dados trimestrais para o período1960q1 até 1990q4, O modelo ARIMA mais simples que inclui uma diferenciaçãoprévia da série original e tanto componentes autoregressivas como de médiasmóveis é a especi�cação ARIMA(1,1,1), Nós podemos ajustar este modelo como comando arima através de:
[]arima wpi, arima(1,1,1)
Podemos ver através dos resultados acima que o coe�ciente AR(1) é 0,874,o coe�ciente MA(1) é -0,412 (sendo que estes dois parâmetros são altamentesigni�cativos) e o desvio padrão do termo de ruido branco εt estimado é 0,725.
2
Este modelo também pode ser estimado usando o comando:arima D.wpi, ar(1) ma(1)Neste último caso estamos previamente aplicando o operador diferença D. a
série para depois modelar a série diferenciada como um processo ARMA(1,1).Modelar a série original como um processo ARIMA(1,1,1) é equivalente a mod-elar a série original diferenciada como um processo ARMA(1,1).
Se yt é a série de preços ao atacado dos EUA então estimamos um modeloARIMA com a seguinte equação:
D.yt = 0, 749 + 0, 874D.yt−1 − 0, 412εt−1 + εtcom σεt = 0, 725Exemplo 2: Após examinar as primeiras diferenças da série WPI, Enders
escolhe um modelo de diferenças aplicado ao logaritimo natural da série originalpara estabilizar a variância da série diferenciada. Os dados brutos e a primeiradiferença do logaritimo da séries pode ser representado gra�camente através docomando:
tsset t, quarterlytwoway (tsline wpi), title(Indice de Preços ao Atacado - EUA) name(gra�co1)twoway (tsline D.ln_wpi), title(Indice de Preços ao Atacado - EUA diferença
dos logs) name(gra�co2) yline(0)graph combine gra�co1 gra�co2, commonscheme xsize(10) ysize(5)
2040
6080
100
120
U.S
. Who
lesa
le P
rice
Inde
x
1960q1 1970q1 1980q1 1990q1quarterly date
Indice de Preços ao Atacado - EUA
-.02
0.0
2.0
4.0
6.0
8D
.ln_w
pi
1960q1 1970q1 1980q1 1990q1quarterly date
Indice de Preços ao Atacado - EUA diferença dos logs
Com base nas autocorrelações, autocorrelações parciais (ver grá�cos abaixo)e os resultados de estimações preliminares, Enders identi�cou ummodelo ARMAsob a séries log-diferenciada.
ac D.ln_wpi, ylabels(-.4(.2).6) name(gra�co3)pac D.ln_wpi, ylabels(-.4(.2).6) name(gra�co4)graph combine gra�co3 gra�co4, commonscheme xsize(10) ysize(5)
3
-0.4
0-0
.20
0.00
0.20
0.40
0.60
Aut
ocor
rela
tions
of D
.ln_w
pi
0 10 20 30 40Lag
Bartlett's formula for MA(q) 95% confidence bands
-0.4
0-0
.20
0.00
0.20
0.40
0.60
Par
tial a
utoc
orre
latio
ns o
f D.ln
_wpi
0 10 20 30 40Lag
95% Confidence bands [se = 1/sqrt(n)]
Em adição ao termo autoregressivo e ao ao termo MA(1), um termo MA(4)é incluido para levar em conta o efeito trimestral remanescente. Portanto, omodelo a ser ajustado é:
∆ln(wpit) = β0 + ρ1{∆ln(wpit−1 − β0) + θ1εt−1 + θ4εt−4 + εtPodemos estimar este modelo com o comando arima:arima D.ln_wpi, ar(1) ma(1 4)
4
3 Modelos ARIMA com sazonalidade multiplica-
tiva
Muitas séries temporais exibem um componente periódico sazonal e o modeloARIMA sazonal comumente chamado de SARIMA pode ser usado. Por exem-plo, dados de vendas mensais de aparelhos de ar condicionado tem uma fortecomponente sazonal, com vendas altas nos meses de verão e baixas nos mesesde inverno. No exemplo anterior levamos em conta efeitos trimestrais ao ajustarum modelo:
(1− ρ1L){∆ln(wpit)− β0} = (1 + θ1L+ θ4L4)εt (2)
Este é um modelo ARIMA com sazonalidade aditiva, no sentido de queos termos MA de primeira e de quarta ordem são especi�cados aditivamente
5
na equação (1 + θ1L + θ4L4). Outra forma de considerar o efeito trimestral
é ajustar um modelo ARIMA com sazonalidade multiplicativa. Um modeloSARIMA multiplicativo de ordem (1, 1, 1)(0, 0, 1)4para a série ln(wpit) é:
(1− ρ1L){∆ln(wpit)− β0} = (1 + θ1L)(1 + θ4,1L4)εt (3)
ou após expandir os termos:
∆ln(wpit) = β0 +ρ1{∆ln(wpit−1)−β0}+ θ1εt−1 + θ4,1εt−4 + θ1θ4εt−5 + εt (4)
Na notação (1, 1, 1) × (0, 0, 1)4, o termo (1, 1, 1) signi�ca que há um termoautoregressivo não sazonal (1− ρ1L) e um termo de médias móveis não sazonal(1 + θ1L) e que a série temporal é diferenciada uma vez no tempo. O termo(0, 0, 1)4 indica que não há um termo autoregressivo sazonal, indica que há umtermo de médias móveis sazonal (1 + θ4,1L
4) e que a série não é diferenciadasazonalmente. Este é conhecido como modelo SARIMA multiplicativo porqueos fatores sazonais e não sazonais trabalham de forma multiplicativa: (1 +θ1L)(1 + θ4,1L
4). Multiplicando-se os termos impomos restrições não linearesaos parâmetros dos termos de defasagem de quinta ordem dos valores da sérietemporal; o comando arima impõe estas restrições automaticamente.
Para melhor esclarecer esta notação, considere um modelo multiplicativoSARIMA (2, 1, 1)× (1, 1, 2)4:
(1− ρ1L− ρ2L2)(1− ρ4,1L4)∆∆4zt = (1 + θ1L)(1 + θ4,1L4 + θ4,2L
8)εt (5)
onde ∆denota o operador diferença ∆yt = yt− yt−1 e ∆s denota o operadorde diferença sazonal de lag s ∆syt = yt − yt−s
Expandindo a última expressão temos:
zt = ρ1zt−1+ρ2zt−2+ρ4,1zt−4−ρ1ρ4,1zt−5−ρ2ρ4,2zt−6+θ1εt−1+θθ4,1εt−4+θ1θ4,1εt−5+θ4,2εt−8+θ1θ4,2εt−9+εt(6)
onde zt = ∆∆4zt = ∆(zt − zt−4) = (zt − zt−1)− (zt−4 − zt−5)e zt = yt − xtβ se regressores são incluidos no modelo.Mais geralmente um modelo multiplicativo SARIMA (p, d, q)× (P,D,Q)sé:ρ(L
p
)ρs(LP )∆d∆D
s zt = θ(Lq)θs(LQ)εt
onde:ρs(L
P ) = (1− ρs,1Ls − ρs,2L2s − ...− ρs,PLPsθ(LQ) = (1 + θs,1L
s + θs,2L2s − ...− θs,Qρs,PLQs
Se o grá�co dos dados sugerem que o efeito sazonal é proporcional a médiada série, então o efeito sazonal é provavelmente multiplicativo e um modeloSARIMA multiplicativo pode ser apropriado. Box, Jenkins, and Reinsel (2008,sec. 9.3.1) sugerem começar com um modelo SARIMA com quaisquer dadosque exibirem um padrão sazonal e depois explorar modelos não multiplicativosSARIMA se os modelos multiplicativos não se ajustarem bem aos dados. Por
6
outro lado, Chat�eld (2004, 14) sugere que tomar logaritimo da série tornará oefeito sazonal aditivo e neste caso um modelo SARIMA aditivo como o ajustadono exemplo anterior seria apropriado. Em resumo, o analista deve tentar tantoo SARIMA aditivo como o SARIMA multiplicativo e veri�car qual dos doisfornece melhor ajuste aos dados e melhores previsões.
Exemplo 3: Modelo SARIMA multiplicativoUma das mais comuns especi�cações SARIMA é o (0, 1, 1)× (0, 1, 1)12 apli-
cado aos dados de linhas aéreas de Box, Jenkins, and Reinsel (2008, sec. 9.2).O conjunto de dados airline.dta contem dados mensais de passageiros de lin-has aéreas internacionais de janeiro de 1949 a dezembro de 1960. Após umaprimeira diferenciação e uma diferenciação sazonal não podemos suspeitar deuma componente de tendência nos dados e então usamos a opção nocontant como comando arima:
use http://www.stata-press.com/data/r13/air2generate lnair = ln(air)arima lnair, arima(0,1,1) sarima(0,1,1,12) noconstant
7
[]Portanto, o nosso modelo estimado de número mensal de passgeiros delinhas internacionais é:
∆∆12lnairt = −0, 402εt−1 − 0, 557εt−12 + 0, 224εt−13 + εtσε = 0, 037O coe�ciente de εt−13(0,224) é aproximadamente igual ao produto dos coe-
�cientes deεt−1(−0, 402) e εt−12(−0, 557) . O comando arima designou o termoDS12.lnair para indicar que foi aplicado o operador diferença ∆e o operadorde diferença sazonal de 12 lags ∆12 a série lnair. Podemos também ajustar omesmo modelo através do comando:
arima DS12.lnair, ma(1) mma(1, 12) noconstantPara modelos multiplicativos SARIMA mais simples é mais fácil usar a opção
sarima(), embora esta segunda sintaxe permita incorporar termos sazonais maiscomplicados.
8
4 Previsões dinâmicas após o arima
Outra característica do comando arima é sua habilidade para utilizar o comandopredict logo após a estimação do modelo para fazer previsões dinâmicas para asérie temporal. Suponha que desejamos ajustar o modelo de regressão:
yt = β0 + β1xt + ρyt−1 + εtusando uma amostra de dados det = 1, .., T e fazer previsões começando no
tempo f . Se usamos o comando regress para ajustar o modelo, então podemosusar o comando predict para fazer previsões um passo a frente. Ou seja, ocomando predict calculará:
yf = β0 + β1xf + ρyf−1
O que é importante destacar aqui é que o comando predict usará o valorefetivo de y no períodof − 1para calcular a previsão no período f . Portanto, seusarmos o comando regress não podemos fazer previsões para períodos adiantede f = T + 1 a menos que tenhamos valores observados para y para estesperíodos.
Se ao invés disso ajustrarmos nosso modelo através do comando arima, entãoo comando predict pode preduzir previsões dinâmicas através do uso do �ltrode Kalman. Se usamos a opção dynamic(f) , então para o período f o comandopredict calculará:
yf = β0 + β1xf + ρyf−1
através do uso do valor observado de yf−1 justamente como quando usamoso comando predict após o comando regress. Entretando, para o período f + 1 ocomando predict < newvar >, dynamic(f) calculará:
yf+1 = β0 + β1xf+1 + ρyfusando o valor predito de yf ao invés do valor observado. De forma similar,
para o período f + 2, a previsão será:yf+2 = β0 + β1xf+2 + ρyf+1
Naturalmente, como o nosso modelo inclui o regressor xt, podemos fazerprevisões apenas para períodos para os quais temos observações para xt. En-tretanto, para modelos ARIMA puros, podemos calcular previsões dinãmicasmuito além do período �nal de nosso conjunto de dados se assim desejarmos.
No exemplo 2 ajustamos o seguinte modelo:∆ln(wpit) = β0 + ρ1{∆ln(wpit−1)− β0}+ θ1εt−1 + θ4εt−4 + εtatravés dos comandos:use http://www.stata-press.com/data/r13/wpi1arima D.ln_wpi, ar(1) ma(1 4)(resultados omitidos)predict xb,xbE o Stata calcula xbtcomo:xbt = β0 + ˆρ1{∆ln(wpit−1)− β0}+ θ1εt−1 + θ4εt−4 + εtonde
εt−j =
{∆ln(wpit−j)− xbt−j t− j > 0
0 em caso contrario
9
signi�cando que o comando predict newvar,xb calcula predições usando amétrica da variável dependente. Neste exemplo, a variável dependente é rep-resentada pelas mudanças em ln(wpit) e portanto, as predições são, da mesmaforma, mudanças nesta variável.
Se, ao invés, usamos:predict y, yStata calcula yt como yt = xbt+ ln(wpit) de forma que ytrepresenta os níveis
preditos de ∆ln(wpit). Em geral predict newvar, y reverterá qualquer operadorde série temporal aplicado a variável dependente durante a estimação.
Se nós queremos ignorar os componentes de erro ARMA quando fazemos aspredições, usamos a opção structural:
predict xbs, xb structuralque gera xbst = β0 porque não há regressores no modelo epredict ys, y structuralgera yst = β0 + ln(wpit−1)Exemplo 4: Previsões dinâmicasUma característica atrativa do comando arima é sua habilidade para fazer
previsões dinâmicas. Suponhamos que queiramos estimar uma função consumocomo:
consumpt = β0 + β1m2t + µtµt = ρµt−1 + θεt−1 + εt
Primeiro, nós estimaremos o modelo usando os dados até o primeiro trimestesde 1978 e então calcularemos as previsões um passo adiante e as previsõesdinâmicas.
use http://www.stata-press.com/data/r13/friedman2keep if time<=tq(1981q4)arima consump m2 if tin(, 1978q1), ar(1) ma(1)(resultados omitidos)Para fazer as previsões um passo adiante:predict chat, yPorque nossa variável dependente não contem operadores de séries temporais
poderíamos ter usado predict chat, xb e teriamos obtido os mesmos resultados.Faremos também previsões dinâmicas, mudando dos valores observados de con-sump para os valores previstos no primeiro trimestre de 1978:
predict chatdy y, dynamic(tq(1978q1)) ytwoway (line consump time, lwidth(medthick) lpattern(dot)) (line chat time,
lpattern(dash)) (line chatdy time) if time>=tq(1977q1), ytitle(Bilhões de dólares)xtitle(trimestre)
10
1200
1400
1600
1800
2000
Bilh
ões
de d
ólar
es
1977q1 1978q1 1979q1 1980q1 1981q1 1982q1trimestre
série observada previsão um passo a frenteprevisão dinâmica
5 Aplicação de modelos ARIMA para séries de
preços
use indices de preços.dta, cleargen time = _ntsset timetwoway (line Vestuario time), xlabel(#44, labsize(vsmall) grid)gen dvestuario = Vestuario - L12.Vestuariotwoway (line dvestuario time), xlabel(#44, labsize(vsmall) grid)
11
5010
015
020
0V
estu
ario
0 2 4 6 8 10 1214 16 1820 22 2426 2830 32 3436 38 4042 4446 48 5052 54 5658 6062 64 6668 70 7274 7678 80 8284 86 88
time
A série apresenta um padrão cíclico com padrão nitidamente sazonal. Éimportante observar que os picos de preço ocorrem sempre no mes de dezembro,sugerindo um processo de formação de preços nitidamente in�uenciado pelademanda. Os pontos mais baixos da série costumam ocorrer no mês de fevereiro.A série (de�acionada) parece também ser estacionária. Podemos veri�car istoatravés do comando <dfuller Vestuario>. Vamos tentar modelar a série usandoinicialmente um modelo SARIMA (0, 0, 1)×(1, 1, 10)12 e depois tentando outrasespeci�cações:
ac dvestuariopac dvestuarioarima Vestuario, arima(0,0,1) sarima(1,1,0,12) noconstantestimates store modelo1arima Vestuario, arima(1,0,0) sarima(1,1,0,12) noconstantestimates store modelo2arima Vestuario, arima(0,0,0) sarima(1,1,0,12) noconstantestimates store modelo3arima Vestuario, arima(1,0,1) sarima(1,1,0,12) noconstantestimates store modelo4estimates stats modelo1 modelo2 modelo3 modelo4estimates restore modelo4ac dvestuariopac dvestuarioarima Vestuario, arima(0,0,1) sarima(1,1,0,12) noconstantestimates store modelo1arima Vestuario, arima(1,0,0) sarima(1,1,0,12) noconstant
12
estimates store modelo2arima Vestuario, arima(0,0,0) sarima(1,1,0,12) noconstantestimates store modelo3arima Vestuario, arima(1,0,1) sarima(1,1,0,12) noconstantestimates store modelo4estimates stats modelo1 modelo2 modelo3 modelo4
Note: N=Obs used in calculating BIC; see [R] BIC note
modelo4 76 . -254.9083 4 517.8166 527.1395
modelo3 76 . -264.7449 2 533.4899 538.1513
modelo2 76 . -258.7182 3 523.4363 530.4285
modelo1 76 . -260.9831 3 527.9662 534.9584
Model Obs ll(null) ll(model) df AIC BIC
O último modelo (modelo4) parece ser o que pelos critérios (máxima verossim-ilhança, AIC e BIC) o que melhor se ajusta aos dados.
estimates restore modelo4predict ruido, residualstwoway (line ruido time)ac ruido
-0.4
0-0
.20
0.00
0.20
0.40
Aut
ocor
rela
tions
of r
uido
0 10 20 30 40Lag
Bartlett's formula for MA(q) 95% confidence bands
predict chat, ypredict chatdy, dynamic(70) ylabel var Vestuario "série observada"label var chat "previsão um passo a frente"label var chatdy "previsão dinâmica"
13
5010
015
020
0ín
dice
de
preç
os v
estu
ário
60 70 80 90mês
série observada previsão um passo a frenteprevisão dinâmica
14
