MODELO LINEAR GENERALIZADO AUTORREGRESSIVO MÉDIA MÓVEL

Exemplo 1: Número de casos de poliomielite (Polio) reportados pelo Centro de Controle de Doenças

dos Estados Unidos, de janeiro de 1970 a dezembro de 1983. Esta é uma série de contagens com valores

baixos, portanto a suposição de distribuição Normal não seria adequada neste caso.
14
12
10
8
Polio

6
4
2
0

1970 1972 1974 1976 1978 1980 1982 1984

Figura 1: Série Polio

Se estamos interessados em modelar séries temporais de contagem, como a série do Exemplo 1,

não podemos utilizar os modelos ARIMA de Box e Jenkins, já que neste caso temos contagens muito

baixas e as observações não têm distribuição normal.

Para resolver este problema, uma possibilidade é a utilização de modelos lineares generalizados

(MLG), pois esta classe de modelos permite o uso de várias distribuições de probabilidade que pertençam

à família exponencial. Por exemplo, séries de contagem podem ter distribuição binomial, Poisson ou

binomial negativa. Entretanto estes modelos, apesar de não necessitarem da suposição de normalidade,
são baseados na suposição de independência das observações, o que geralmente não ocorre em séries

temporais.

Assim, surgiram os modelos lineares generalizados com erros autorregressivos e médias móveis

(GLARMA), propostos por Davis et al. (2003). Estes modelos utilizam o MLG conjuntamente com o

ARMA para modelar uma gama muito maior de séries temporais, que tenham qualquer distribuição

pertencente à família exponencial.

1. Modelos Lineares Generalizados

Os modelos lineares generalizados (MLG) são uma extensão dos modelos lineares normais e

foram propostos por Nelder e Wedderburn (1972). A ideia básica consiste em abrir o leque de opções

para a distribuição da variável resposta, permitindo que a mesma pertença à família exponencial de

distribuições.

Seja y1, . . . , yn uma amostra aleatória de n observações independentes de uma variável resposta

Y. No modelo linear geral, supomos que Yt tem distribuição normal com média µt e variância σ2,

( )
Yt ~ N t , 2 . Além disto, assumimos que o valor esperado, µt, é uma função linear de r variáveis

preditoras x ~ (x1 ,..., xr ) , i.e.,

E (Yt ) =  t = x't  ,

onde β é um vetor de parâmetros desconhecidos.

Mas vamos supor que a distribuição de Y não seja Normal. Precisamos determinar a distribuição

correta da variável resposta para podermos estimar os parâmetros, fazer testes de hipóteses e tirar

conclusões sobre o modelo. Uma classe de distribuições que garante essas propriedades é conhecida

como família exponencial.

2
 1.1. Família exponencial a um parâmetro

A família exponencial (FE) é uma família de distribuições cuja função densidade pode ser escrita

na seguinte forma:

f ( yt ,t ) = expa( yt )b(t ) + c(t ) + d ( yt ), t = 1,..., n (1.1)

onde t é o parâmetro de interesse e a(.), b(.), c(.) e d(.) são funções conhecidas.

Obs.: Se a(y) = y então b(.) é chamado parâmetro canônico da distribuição.

Se Yt tem distribuição na família exponencial, então sua média e variância são dadas por

c' (t )
E (Yt ) = t = −
b' (t )

b' ' (t )c' (t ) − c' ' (t )b' (t )
Var (Yt ) = .
b' (t )3

1.2. Componentes do Modelo

Uma vez definido como se estruturam as distribuições segundo a família exponencial, temos três

componentes que formam um MLG:

a) Componente Aleatória – composto de uma variável resposta Y com distribuição pertencente

à família exponencial;

b) Componente Sistemática – é a parte do modelo composta pelas variáveis explicativas

x = ( x1 ,..., xr ) e denotada por  t = x t'  , onde ηt é conhecido como preditor linear.

c) Função de Ligação (g) – é uma função diferenciável e contínua. Através desta função, faz-se

a conexão entre a média das observações e a parte sistemática:

3
 g (t ) =  t = x t' 

Alguns exemplos de funções de ligação são dados na Tabela 1.

Tabela 1: Alguns exemplos de funções de ligação

Distribuição Ligação

Normal Identidade

Poisson Log

Binomial Logit

1.3. Estimação

O vetor de parâmetros β = (β1, ..., βr) pode ser estimado através de várias metodologias, como

mínimos quadrados ponderados, máxima verossimilhança (MV) ou inferência Bayesiana. Para o método

MV, a função de log-verossimilhança é dada por:

n
( , y ) =  (a( yt )b(t ) + c(t ) + d ( yt ) ) . (1.2)
t =1

Como as derivadas da função acima são difíceis de serem obtidas analiticamente, devemos usar

métodos numéricos para encontrar as estimativas dos parâmetros. Um dos procedimentos mais utilizados

é o processo de otimização iterativo conhecido como algoritmo de Newton-Raphson. Resumidamente,

pode-se dizer que o algoritmo inicia o processo especificando uma estimativa inicial para β e vai

sucessivamente alterando-a até que a diferença entre β na iteração (m + 1) e a estimativa anterior seja

menor que um valor bem pequeno pré-definido. Desta forma, a convergência é obtida e temos o vetor de

coeficientes estimados.

4
 Através das propriedades do estimador de máxima verossimilhança, para n suficientemente

grande temos,

( )
ˆ ~ N r  , I −1 ,

onde I é a matriz de informação de Fisher.

1.4. Adequação de modelos

Após a estimação dos coeficientes do modelo, alguns procedimentos devem ser realizados a fim

de medir a qualidade do ajuste e adequabilidade do modelo, ou seja, deve ser feita uma validação dos

resultados. Num primeiro momento essas análises servem também para comparar diferentes modelos.

Assim como no modelo linear normal, para se decidir entre um ou outro modelo, pode-se calcular

medidas associadas ao ajuste do modelo, como o coeficiente de determinação (R²) e o critério de

informação AIC.

Além destas medidas, uma outra maneira de verificar a adequação de um modelo é compará-lo

com um modelo mais geral, com o número máximo de parâmetros que podem ser estimados. Este último

é chamado de modelo saturado. A estatística Desvio (do inglês Deviance) calcula a bondade de ajuste do

modelo através das diferenças entre a função de log-verossimilhança do modelo saturado com o modelo

sob investigação, i.e., (ˆsat , y ) − (ˆ, y ) . Aqui ̂sat denota o EMV do vetor de parâmetros do modelo

saturado e ˆ é o EMV dos parâmetros do modelo sob investigação.

Assim, a estatística Desvio é definida como:


D = 2 (ˆsat , y ) − (ˆ, y ) . (1.3)

5
 Um valor pequeno desta estatística implica um bom ajuste. A distribuição amostral de D é

aproximadamente  2 (n − r , ) , onde  é o parâmetro de não-centralidade.

1.5. Análise de resíduos

No caso do MLG, existem dois tipos de resíduos que são mais utilizados na prática, os resíduos

de Pearson e Desvio.

a) Resíduo de Pearson

O resíduo de Pearson é simplesmente a diferença entre os valores observados e ajustados, dividido

por uma estimativa do desvio padrão do valor ajustado. O resíduo resultante tem a forma,

yt − ˆ t
rtP = .
Var (ˆ t )

b) Resíduo Desvio

Se a estatística Desvio dada na Equação (1.3) é usada como uma medida de discrepância, então

cada unidade contribui com uma quantidade para o Desvio, logo D =  d t . Desta forma, o resíduo
t

Desvio é definido como:

rtD = sinal( yt − t ) dt .

Ou seja, se tomarmos o quadrado de rtD e somarmos sobre todas as observações, obtemos a

estatística Desvio. Observações com um resíduo Desvio maior que 2 podem indicar falta de ajuste.

6
 Pode-se construir gráficos dos resíduos versus tempo e observar se eles se encontram

aleatoriamente distribuídos em torno de zero, com variância constante. Além disto, deve-se também

observar se há indícios de autocorrelação nos resíduos, construindo o gráfico da FAC.

1.6. Exemplo – Série Polio

O Exemplo 1 mostrou a série Polio, referente ao número de casos de poliomielite nos Estados

Unidos, de janeiro de 1970 a dezembro de 1983. A Figura 2 apresenta o histograma dos dados. Podemos

verificar que não é possível utilizar a distribuição Gaussiana neste caso.

Histogram of polio
140
120
100
80
Frequency

60
40
20
0

0 2 4 6 8 10 12 14

polio

Figura 2. Histograma da série Polio

Muitos autores já analisaram esta série na literatura para verificar se a incidência de poliomielite

vem decrescendo desde 1970. Uma possibilidade para a distribuição da variável resposta, como se trata

de contagens, é a distribuição de Poisson. Como variáveis explicativas, a maioria dos estudos utiliza um

componente de tendência e componentes sazonais usando pares de seno e cosseno, com ciclos anuais e

7
semianuais. Retirando as últimas 12 observações (jan/83 a dez/83) para fazer previsões, o modelo é dado

por:

yt ~ Poisson (  t ), t = 1,...,156

 2 t   2 t   2 t   2 t 
t = log(t ) =  0 + 1t +  2 cos  +  3 sen  +  4 cos  +  5 sen 
 12   12   6   6 

Saída do R:

Estimate Std. Error z value Pr(>|z|)

(Intercept) 0.622557 0.129897 4.793 1.65e-06 ***

Ano -0.005820 0.001571 -3.704 0.000212 ***

Cos12 0.095782 0.093143 1.028 0.303790

Sen12 -0.517376 0.116958 -4.424 9.71e-06 ***

Cos6 0.405836 0.103904 3.906 9.39e-05 ***

Sen6 -0.088790 0.101068 -0.879 0.379661

Null deviance: 318.76 on 155 degrees of freedom

Residual deviance: 271.09 on 150 degrees of freedom

AIC: 522.87

O coeficiente de 1 (Ano) foi negativo e significativo, o que significa que a incidência de pólio

está diminuindo com o tempo. Além disto, vemos que os componentes sazonais de sen12 e cos6 são

significativos, portanto devemos manter também os componentes cos12 e sen6.

Assim, o modelo ajustado é dado por:

  2 t   2 t   2 t   2 t  
ˆ t = exp 0,623 − 0,006 t + 0,096 cos  − 0,517 sen  + 0,406 cos  − 0,089sen  
  12   12   6   6 

8
 A Figura 3 apresenta os gráficos de resíduos versus tempo e a FAC, tanto para o resíduo de

Pearson, quanto o Desvio. O gráfico de resíduos no tempo parece apresentar um comportamento aleatório

em torno do valor zero, apesar de termos algumas observações com valores muito altos (acima de 3

desvios-padrão). A FAC também mostra a falta de independência dos resíduos, pois vemos alguns picos

significativos nos lags iniciais.

Resíduos de Pearson Series rp

0.2 0.4 0.6 0.8 1.0

6
4

ACF
rp

2
0

-0.2

0 50 100 150 0 5 10 15 20

Index Lag

Resíduos Desvio Series rd

4 5

0.2 0.4 0.6 0.8 1.0

-2 -1 0 1 2 3

ACF
rd

-0.2

0 50 100 150 0 5 10 15 20

Index Lag

Figura 3. Análise de resíduos para o modelo [Link]

9
Previsões: Vamos fazer previsões 12 passos à frente utilizando este modelo.

Janeiro de 1983: A previsão da incidência de poliomielite para o mês de jan/83 é:

  314   314   314   314 

exp 0,623 − 0,006 x157 + 0,096 cos  − 0,517 sen  + 0,406 cos  − 0,089sen   = 0,71
  12   12   6   6 

Desta forma, podemos calcular as previsões para o período Jan/83 a Dez/83, que são apresentadas

na Tabela 2. A Figura 4 mostra o ajuste, assim como previsões para os doze últimos meses. Podemos ver

que tanto o ajuste quanto as previsões parecem satisfatórios, somente apresentando valores mais

afastados dos reais nos períodos de picos na incidência de pólio.

Tabela 2: Previsões 12 passos à frente para a série de Polio, Jan/83 a Dez/83

Mês Real Previsto Real-Previsto

Jan/1983 0 0,71 -0,71
Fev/1983 1 0,38 0,62
Mar/1983 0 0,29 -0,29
Abr/1983 0 0,39 -0,39
Mai/1983 0 0,69 -0,69
Jun/1983 1 0,99 0,01
Jul/1983 2 0,97 1,02
Ago/1983 1 0,81 0,19
Set/1983 0 0,80 -0,80
Out/1983 1 1,03 -0,03
Nov/1983 3 1,31 1,69
Dez/1983 6 1,16 4,84

EQMP = 2,4687.

10
 10
8
6
Polio

4
2
0
-2

0 50 100 150

tempo

Figura 4: Ajuste e previsão para o modelo [Link]. A linha preta representa a série Polio, a linha azul
mostra o modelo ajustado e os pontos em azul são os valores previstos.

11
2. Modelos GLARMA

Neste capítulo, descrevemos o modelo linear generalizado autoregressivo média móvel

(GLARMA), introduzido por Davis et al. (2003), que é útil para modelar variáveis respostas dependentes

no tempo e que seguem uma distribuição da família exponencial. Daremos especial atenção ao modelo

de Poisson, já que esta é a distribuição mais utilizada para dados de contagem.

2.1. Definição

A classe GLARMA é uma classe de modelos que estende o processo ARMA Gausssiano de séries

temporais para um modelo mais flexível para séries de contagem não-Gaussianas. A variável dependente

é suposta ter uma distribuição condicional na família exponencial dado todo o passado do processo.

( )
Sejam yt a série temporal e Ft −1 = y (t −1) , x (t ) , onde y (t −1) = ( y1 , y2 ,..., yt −1 ) é o passado do

processo de contagem e x (t ) é o passado e o presente das variáveis regressoras.

Para introduzir o modelo GLARMA, assuma que yt, dado o passado histórico Ft −1 , tem qualquer

distribuição na família exponencial (FE) com média  t ,

Yt | Ft −1 ~ FE ( t ) , t = 1,..., n . (2.1)

O preditor  t é dado por:


t = g ( t ) =  0 + 1 x1, t + ... +  r xr , t +  i  t − i (2.2)
i =1

onde  é um vetor de parâmetros e  t é um ruído branco.

12
 
O componente  i  t − i pode ser especificado em termos de um número finito de parâmetros
i =1

utilizando a metodologia de Box e Jenkins (1976),

  q ( B)
 ( B ) =  i B i = −1
i =1  p ( B)

onde  p (B) e  q (B) são, respectivamente, os polinômios autorregressivo e média móvel.

Dessa forma o modelo GLARMA(p,q) fica especificado por

Yt | Ft −1 ~ FE ( t ) , t = 1,..., n

onde t = g (t ) =  0 + 1x1,t + ... +  r xr ,t + Tt =  0 + xt'  + Z t

( )
Z t = 1 (Z t −1 +  t −1 ) + ... +  p Z t − p +  t − p + 1 t −1 + ... +  q t − q e

yt −  t
t = .
(t )
Exemplo:

GLARMA(0,1) com distribuição de Poisson:

Yt | Ft −1 ~ Poisson (  t ) , t = 1,..., n

onde  t = g (  t ) =  0 + x t'  + Z t ,

Z t =  t −1 ,

yt −  t
t = .
(t )

13
Assim,

t =  0 + x t'  +  t −1 =  0 + x t'  + 
( yt −1 − t −1 ) =  + x '  +  ( yt −1 − exp(t −1 ))
(t −1 ) (exp(t −1 ))− 
0 t

Neste caso, E ( t ) = E (E ( t |  t −1 )) =  0 + x t'  .

2.2. Estimação

A estimação dos parâmetros do GLARMA,  = (  ,  )' , onde  = (1 ,..., p ,1 ,..., q ) e

 = (  0 , 1 ,...,  r ) , é feita conjuntamente através da função de verossimilhança, maximizada pelo

método numérico Newton-Raphson (Davis et al., 2003).

O valor  de para o cálculo dos resíduos deve ser especificado pelo pesquisador. Se  = 0,5

temos o resíduo de Pearson. Outro valor usado na prática é  = 1 .

Considere a densidade condicional de Yt dado Ft −1 pertencente à família exponencial. A função

de log-verossimilhança pode ser escrita como

n
( , y ) =  log f ( yt | Ft −1 ) . (2.3)
t =1

No caso da distribuição de Poisson, a função de log-verossimilhança é dada por:

n
( , y ) =  (Ytt − exp(t )) . (2.4)
t =1

14
 Para inicializar o método recursivo de Newton-Raphson na maximização numérica da log-

verossimilhança ( , y) , Davis et al. (2003) sugerem que os valores obtidos das estimativas do

GLARMA sem os termos auto-regressivos média móveis sejam utilizados como valores iniciais. A

convergência, na maioria dos casos, ocorre após 10 iterações. A matriz de covariância dos estimadores é

estimada por

−1
  2 (ˆ , y ) 
 = −
ˆ  .
  ' 
 

Maiores detalhes sobre as condições de estacionariedade, propriedades, estimação e inferência

dos modelos GLARMA podem ser vistos em Davis et al. (2003).

A análise de resíduos no modelo GLARMA é feita da mesma forma que nos modelos GLM e

ARMA. Usando os resíduos estimados

yt − ˆ t
et =
(ˆ t )
devemos verificar se os mesmos estão aleatoriamente distribuídos em torno de zero, se têm variância

constante e se são independentes.

2.3. Previsão

Adaptando a metodologia de previsões de modelos ARMA para o modelo GLARMA, tem-se que

o valor previsto de origem t e horizonte h, YˆT (h ) é dado por:

( )
YˆT (h) = E (YT + h | YT ) = ET (YT + h ) = ˆ T + h = g −1 ˆ0 + x t' ˆ + ET (Z T + h ) (2.5)

15
onde

( ( ) ( )) ( (
ET (Z T + h ) = ˆ1 (ET (Z T + h−1 ) + ET ( T + h−1 )) + ... + ˆ p ET Z T + h− p + ET  T + h− p + ˆ1 (ET ( T + h−1 )) + ... + ˆq ET  T + h−q ))

Para calcular as esperanças condicionais na prática, deve-se considerar:

i) Et Z t − h  = Z t − h para h = 0, 1, 2, ...

ii) Et Z t + h  = Zˆ t (h ) para h = 1, 2, ...

iii) Et  t − h  =  t − h para h = 0, 1, 2, ...

iv) Et  t + h  = 0 para h = 1, 2, ...

Por exemplo, para um GLARMA(1,1) com distribuição de Poisson:

(
YˆT (h) = exp ˆ0 + x t' ˆ + ˆ1 (Et (Z t + h −1 ) + Et ( t + h −1 )) + ˆ1 (Et ( t + h −1 )) .)
Para h = 1:

(
YˆT (1) = exp ˆ0 + x t' ˆ + ˆ1 (Z t +  t ) + ˆ1   t )
Para h ≥ 2:

(
YˆT (h) = exp ˆ0 + x t' ˆ + ˆ1  Zˆ t (h) . )

2.4. Exemplo – Série Polio

Vamos voltar ao exemplo da série Polio. Se supomos que os dados seguem uma distribuição de

Poisson, podemos ajustar o modelo GLARMA-Poisson. Vamos usar o pacote “glarma” do R.

16
 Na Seção 1 vimos, pela FAC e FACP da série Polio, que um possível modelo para esta série é

um AR(1). Desta forma, inicialmente consideramos um modelo GLARMA(1,0).

Modelo GLARMA utilizando λ = 0,5:

M05.AR1: GLARMA(1,0):

Estimate [Link] z-ratio Pr(>|z|)

phi_1 0.21659 0.04688 4.62 3.83e-06 ***

Linear Model Coefficients:

Estimate [Link] z-ratio Pr(>|z|)
Inter 0.538753 0.179353 3.004 0.002666 **
Ano -0.005508 0.002119 -2.599 0.009350 **
Cos12 0.120156 0.120085 1.001 0.317024
Sen12 -0.521244 0.147593 -3.532 0.000413 ***
Cos6 0.439291 0.117429 3.741 0.000183 ***
Sen6 -0.005650 0.114400 -0.049 0.960609

Null deviance: 318.76 on 155 degrees of freedom

Residual deviance: 246.82 on 149 degrees of freedom
AIC: 506.2899

Number of Fisher Scoring iterations: 30

LRT and Wald Test:

Alternative hypothesis: model is a GLARMA process
Null hypothesis: model is a GLM with the same regression structure
Statistic p-value
LR Test 18.58 1.63e-05 ***
Wald Test 21.35 3.83e-06 ***

Todas as variáveis, exceto cosseno anual e seno semianual foram significativas, porém vamos

deixá-las no modelo para modelar a sazonalidade. A saída do pacote “glarma” também mostra os testes

de Wald e Razão de Verossimilhança, que indicam que o modelo MP.AR1 é superior ao GLM.

A sobrefixação deste modelo mostra que não é necessário inserir mais nenhum parâmetro

autorregressivo no modelo. A utilização de um modelo MA(1) também retorna um ajuste pior, com AIC=

508,926.

17
Modelo GLARMA utilizando λ = 1:

M1.AR1: GLARMA(1,0):

Estimate [Link] z-ratio Pr(>|z|)

phi_1 0.28701 0.04815 5.961 2.51e-09 ***

Linear Model Coefficients:

Estimate [Link] z-ratio Pr(>|z|)
Inter 0.500882 0.189894 2.638 0.008347 **
Ano -0.004745 0.002242 -2.116 0.034317 *
Cos12 0.115437 0.124457 0.928 0.353655
Sen12 -0.476728 0.146218 -3.260 0.001113 **
Cos6 0.420584 0.110813 3.795 0.000147 ***
Sen6 0.130264 0.109267 1.192 0.233200

Null deviance: 318.76 on 155 degrees of freedom

Residual deviance: 211.29 on 149 degrees of freedom
AIC: 502.6204

Number of Fisher Scoring iterations: 30

LRT and Wald Test:

Alternative hypothesis: model is a GLARMA process
Null hypothesis: model is a GLM with the same regression structure
Statistic p-value
LR Test 22.25 2.39e-06 ***
Wald Test 35.53 2.51e-09 ***

Utilizando λ = 1, as mesmas variáveis foram significativas, mas este modelo apresenta menor

AIC. A sobrefixação deste modelo mostra que é necessário inserir mais 3 parâmetros autorregressivos,

obtendo assim um AR(4). Este é o modelo que apresenta melhor ajuste, com o menor AIC (486,1959),

porém é o menos parcimonioso, com 3 parâmetros a mais em relação ao M05.AR1. A saída é apresentada

abaixo:

M1.AR4: GLARMA(1,0):

Estimate [Link] z-ratio Pr(>|z|)

phi_1 0.23978 0.04413 5.433 5.54e-08 ***
phi_2 0.14716 0.03624 4.061 4.89e-05 ***
phi_3 -0.16852 0.05509 -3.059 0.00222 **
phi_4 0.09809 0.04405 2.227 0.02595 *

18
 Linear Model Coefficients:
Estimate [Link] z-ratio Pr(>|z|)
Inter 0.509763 0.204007 2.499 0.01246 *
Ano -0.006765 0.002458 -2.752 0.00592 **
Cos12 0.250976 0.129172 1.943 0.05202 .
Sen12 -0.620193 0.129112 -4.804 1.56e-06 ***
Cos6 0.420130 0.106940 3.929 8.54e-05 ***
Sen6 0.034647 0.105405 0.329 0.74238

Null deviance: 318.76 on 155 degrees of freedom

Residual deviance: 307.45 on 146 degrees of freedom
AIC: 486.1959

Number of Fisher Scoring iterations: 30

LRT and Wald Test:

Alternative hypothesis: model is a GLARMA process
Null hypothesis: model is a GLM with the same regression structure
Statistic p-value
LR Test 44.67 4.64e-09 ***
Wald Test 47.10 1.45e-09 ***

A Figura 5 apresenta a análise de resíduos para os modelos M05.AR1 (λ = 0,5) e M1.AR4 (λ =

1). Os resíduos apresentam um comportamento aleatório em torno de zero, com alguns outliers,

principalmente no modelo M1.AR4, mas podemos considerar a variância constante. Com relação à FAC

e FACP, observamos apenas um pico significativo no lag 8 para o modelo M1.AR4, mas isto pode ser

considerado um ruído. Assim, podemos supor que os resíduos são independentes.

19
 0 1 2 3 4 5 Modelo M05.AR1 Modelo M1.AR4

8
6
rp

rs

4
2
0
-2

0 50 100 150 0 50 100 150

Index Index

Modelo M05.AR1 Modelo M1.AR4

1.0

1.0
0.6

0.6
ACF

ACF
0.2

0.2
-0.2

-0.2

0 5 10 15 20 0 5 10 15 20

Lag Lag

Modelo M05.AR1 Modelo M1.AR4

0.00 0.10

0.00 0.10
Partial ACF

Partial ACF
-0.15

-0.15

5 10 15 20 5 10 15 20

Lag Lag

Figura 5: FAC e FACP para os resíduos Pearson e Escore do modelos M05.AR1 e M1.AR4,

respectivamente

20
Previsões:

A Tabela 3 mostra o EQMP para previsões 12 passos à frente utilizando os modelos M05.AR1 e

M1.AR4. Como o EQMP do modelo com λ = 0,5 é menor, vamos apresentar as previsões somente para

este modelo.

Tabela 3: EQMP para os modelos MP2 e MS1

Pearson (M05.AR1) Escore (M1.AR4)

2,4639 2,664557

A Tabela 4 mostra as previsões 12 passos à frente com o modelo M05.AR1 e a Figura 6 mostra

o ajuste e as previsões. Comparando com o ajuste do MLG realizado na Seção 1.6, vemos que o EQMP

do modelo GLARMA é um pouco menor (EQMP_GLARMA=2,4639 e EQMP_GLM=2,4687).

Tabela 4: Previsões 12 passos à frente com o modelo M05.AR1 para a série Polio, Jan/83 a Dez/83

Mês Real Previsto Real-Previsto

Jan/1983 0 0,88 -0,88
Fev/1983 1 0,40 0,60
Mar/1983 0 0,27 -0,27
Abr/1983 0 0,34 -0,34
Mai/1983 0 0,61 -0,61
Jun/1983 1 0,97 0,03
Jul/1983 2 1,01 0,99
Ago/1983 1 0,82 0,18
Set/1983 0 0,75 -0,75
Out/1983 1 0,92 0,08
Nov/1983 3 1,23 1,77
Dez/1983 6 1,19 4,81

21
 14
12
10
8
Polio

6
4
2
0

0 50 100 150

tempo

Figura 6: Ajuste e previsão para o modelo M05.AR1. A linha preta representa a série Polio, a linha azul
mostra o modelo ajustado e os pontos em azul são os valores previstos.

Referências

1. Benjamin, M. A., R. A. Rigby, and D. M. Stasinopoulos (2003). Generalized autoregressive

moving average models. Journal of the American Statistical Association, 98, pp. 214-223.

2. Box, G.E.P. and Jenkins, G.M. (1976) Time Series Analysis: Forecasting and Control. San

Francisco: Holden-Day.

3. Davis, R.A., Dunsmuir, W.T.M. and Streett, S.B. (2003). Observation-driven models for Poisson

counts. Biometrika, 90, pp 777-790.

4. McCullagh, P.; Nelder, J. A. (1989). Generalized Linear Models. Chapman and Hall, London,

second edition.

22
5. Morettin, P.A., Toloi, C.M.C. (2004) Análise de Séries Temporais. São Paulo: Edgard Blucher.

6. Nelder, J, Wedderburn, R. (1972). Generalized Linear Models. Journal of the Royal Statistical

Society. Series A, 135 (3): 370–384.

7. Shumway, R. H. and Stoffer, D. S. (2011) Time Series Analysis and Its Applications: With R

Examples. Springer, New York, third edition

23
 ANEXO: Séries utilizadas nos exemplos.

Série 1 (TempMedia): Temperatura global de 1900-1997 (calculada como desvio da temperatura global

média anual do período 1961-1990).

-0.13 -0.22 -0.37 -0.44 -0.49 -0.37 -0.30 -0.50 -0.52 -0.49 -0.46 -0.49

-0.41 -0.42 -0.24 -0.13 -0.36 -0.51 -0.39 -0.30 -0.23 -0.19 -0.30 -0.27

-0.33 -0.22 -0.08 -0.19 -0.22 -0.37 -0.13 -0.05 -0.10 -0.23 -0.11 -0.15

-0.10 0.00 0.10 0.02 -0.04 0.06 0.06 0.06 0.22 0.06 -0.08 -0.08

-0.08 -0.09 -0.19 -0.05 0.02 0.10 -0.15 -0.16 -0.26 0.05 0.12 0.04

0.00 0.03 0.04 0.07 -0.22 -0.16 -0.06 -0.06 -0.09 0.03 -0.03 -0.19

-0.06 0.08 -0.18 -0.12 -0.22 0.06 -0.03 0.06 0.10 0.14 0.05 0.24

0.02 0.00 0.09 0.23 0.25 0.18 0.35 0.29 0.15 0.19 0.26 0.39

0.22 0.43

Série 2 (FARELO): Preço do farelo de soja, no estado de São Paulo, no período de jan/1990 a set/1999

Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

1990 278.6 225.1 205.7 161.5 179.1 161.9 163.0 164.8 171.6 175.2 188.4 187.5

1991 187.4 180.0 159.1 153.0 150.1 154.6 158.7 185.2 191.3 208.6 188.6 204.3

1992 241.7 204.4 177.4 159.1 154.7 168.4 183.0 198.4 227.6 210.5 196.5 206.1

1993 209.5 174.8 151.5 149.2 154.3 179.9 228.7 221.0 207.3 207.7 216.3 201.8

1994 192.2 174.2 167.4 153.9 153.2 171.1 172.1 175.0 175.9 181.5 181.5 174.0

1995 178.8 183.4 135.1 127.9 130.6 133.7 158.8 162.9 167.4 181.4 208.1 221.6

1996 240.6 215.1 196.0 219.4 227.3 223.8 234.4 248.3 298.4 275.1 290.1 281.7

1997 277.0 276.6 283.4 282.8 282.1 262.1 243.9 273.6 290.8 270.0 267.6 259.4

1998 239.1 201.1 159.6 153.1 150.7 148.8 147.5 153.9 151.3 154.8 158.3 157.5

1999 152.0 158.1 129.4 130.2 126.3 131.2 133.8 144.8

Série 3 (GRAO): Preço do grão de soja, no estado de São Paulo, no período de jan/1990 a set/1999

24
 Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

1990 12.63 11.04 10.86 9.77 10.55 9.85 9.94 10.10 9.86 9.90 10.16 9.32

1991 9.19 10.02 10.03 9.95 9.85 9.62 9.12 10.10 11.10 11.98 10.05 9.53

1992 11.24 10.55 10.19 9.41 9.72 10.28 10.22 11.02 12.45 11.74 11.86 11.79

1993 12.37 11.07 10.07 9.65 9.80 9.86 12.54 12.81 11.77 11.89 12.71 12.93

1994 12.77 12.15 11.52 10.58 10.91 11.63 11.49 11.52 11.90 12.39 12.62 12.50

1995 12.77 12.14 9.43 9.92 9.37 9.20 9.69 10.83 10.71 11.25 12.46 13.02

1996 15.08 13.08 12.07 13.17 13.86 13.05 13.47 14.36 16.76 16.37 16.57 16.58

1997 15.07 14.56 14.97 15.58 15.98 14.98 14.61 15.70 16.93 16.58 17.37 17.40

1998 15.48 13.02 11.63 11.07 11.45 11.35 11.31 10.98 11.44 11.55 11.43 11.00

1999 9.60 8.65 8.32 8.51 8.68 8.72 8.49 9.09

Série 4 (CEP): Consumo de energia elétrica das Centrais Elétricas do Paraná (CEP), de jan/80 a dez/92

Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

1980 256 261 275 283 293 290 281 292 289 291 296 286

1981 270 283 285 297 303 302 298 301 302 301 290 288

1982 279 295 293 306 313 314 307 311 313 310 311 313

1983 289 305 318 325 332 324 322 334 335 330 320 306

1984 284 304 331 351 365 353 352 354 346 343 321 318

1985 304 337 343 358 364 363 357 361 358 359 329 337

1986 314 356 357 371 383 375 367 368 378 372 338 340

1987 316 361 366 388 395 403 391 394 403 389 369 365

1988 345 383 400 406 428 424 422 426 423 420 392 396

1989 373 407 413 430 443 446 444 450 448 447 417 411

1990 387 422 429 444 450 451 456 455 452 443 420 423

1991 408 438 464 470 478 482 469 471 474 476 452 451

1992 425 465 474 485 506 499 481 492 514 515 483 481

25
Série 5 (Acidente) Número total de motoristas mortos ou feridos na Grã Bretanha devido a acidentes

de trânsito entre Jan/1969 a Dez/1984

Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

1969 1687 1508 1507 1385 1632 1511 1559 1630 1579 1653 2152 2148

1970 1752 1765 1717 1558 1575 1520 1805 1800 1719 2008 2242 2478

1971 2030 1655 1693 1623 1805 1746 1795 1926 1619 1992 2233 2192

1972 2080 1768 1835 1569 1976 1853 1965 1689 1778 1976 2397 2654

1973 2097 1963 1677 1941 2003 1813 2012 1912 2084 2080 2118 2150

1974 1608 1503 1548 1382 1731 1798 1779 1887 2004 2077 2092 2051

1975 1577 1356 1652 1382 1519 1421 1442 1543 1656 1561 1905 2199

1976 1473 1655 1407 1395 1530 1309 1526 1327 1627 1748 1958 2274

1977 1648 1401 1411 1403 1394 1520 1528 1643 1515 1685 2000 2215

1978 1956 1462 1563 1459 1446 1622 1657 1638 1643 1683 2050 2262

1979 1813 1445 1762 1461 1556 1431 1427 1554 1645 1653 2016 2207

1980 1665 1361 1506 1360 1453 1522 1460 1552 1548 1827 1737 1941

1981 1474 1458 1542 1404 1522 1385 1641 1510 1681 1938 1868 1726

1982 1456 1445 1456 1365 1487 1558 1488 1684 1594 1850 1998 2079

1983 1494 1057 1218 1168 1236 1076 1174 1139 1427 1487 1483 1513

1984 1357 1165 1282 1110 1297 1185 1222 1284 1444 1575 1737 1763

26