Trabalho sobre series temporais

Prévia do material em texto

Atividade 5.
2. Carregando o pacote ’fpp2’ no R nos permite utilizar a série ’ausbeer’. Segundo a própria
descrição do R, essa série se trata da produção total de cerveja trimestral na Austrália (em
megalitros) do primeiro trimestre de 1956 até o segundo trimestre de 2010. Considerando
apenas os dados de todos os trimestres de 1960 até 1975, faça o que se pede:
i. Faça o gráfico da série e da sua ACF (editando os eixos). Aparentemente, existe tendência
ou sazonalidade?
library(fpp2)
data()
ausbeer
aust =ts(ausbeer[17:80],frequency=4,start=1960)
aust
autoplot(aust,xlab = ’Tempo (trimestres)’,
ylab = ’produç~ao total de cerveja ’)
acf(aust, main = ’FAC da serie’, xlab = ’Defasagem’,
ylab = ’FAC’)
Figura 1: produção total de cerveja trimestral na Austrália (em megalitros) de 1960 até 1975.
(a) Série temporal (b) FAC da série
Aparentemente, há ind́ıcios da presença de sazonalidade e tendência crescente.
ii. Teste a existência de tendência na série pelo teste de corridas. O que se pôde concluir?
library(DescTools)
RunsTest(aust)
Runs Test for Randomness
data: aust
z = -4.1581, runs = 16, m = 32, n = 32, p-value
= 3.209e-05
alternative hypothesis: true number of runs is not equal the expected number
sample estimates:
median(x)
362
1
Como p-valor< 0, 05, rejeitamos H0 . Dessa forma, há ind́ıcios para acreditar que existe
tendência.
iii. Verifique se existe sazonalidade pelo teste de Kruskal-Wallis (Dica: no script fornecido em
aula, os grupos estão de acordo com os meses, neste caso os trimestres deverão ser utilizados
como grupos). Qual a conclusão?
trimestre = rep(c(’1trim’, ’2trim’, ’3trim’, ’4trim’),16)
length(aust)
[1] 64
length(trimestre)
[1] 64
kw_db = cbind(aust, trimestre)
kruskal.test(aust ~ trimestre, data = kw_db)
Kruskal-Wallis rank sum test
data: aust by trimestre
Kruskal-Wallis chi-squared = 15.259, df = 3,
p-value = 0.001608
Como p-valor< 0, 05, rejeitamos H0 . Dessa forma, há ind́ıcios para acreditar que existe sazo-
nalidade.
iv. Com base nos dois testes, utilize o método de suavização exponencial apropriado.
Identifique o REQM, AIC, BIC e EAM. Faça o gráfico contendo a série original, a série
suavizada e a sua projeção para as próximas 4 observações (4 trimestres de 1976).
fit1 <- hw(aust, seasonal = "additive", h = 4)
fit1
names(fit1)
fit1$model
summary(fit1)
fit1$fitted
plot.ts(aust)
autoplot(aust) +
autolayer(fit1$fitted, series = ’Ajustados’) +
autolayer(fit1, series = ’Projeç~ao’)
2
Utilizando o modelo aditivo de Holt W.
summary(fit1)
Forecast method: Holt-Winters’ additive method
Model Information:
Holt-Winters’ additive method
AIC AICc BIC
608.7296 612.0629 628.1596
Error measures:
ME RMSE MAE MPE MAPE MASE ACF1
Training set -0.1071805 12.62443 9.374296 -0.3158251 2.538721 0.5816523 -0.1164241
Forecasts:
Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
1976 Q1 499.2172 481.9213 516.5132 472.7654 525.6691
1976 Q2 454.7277 437.0609 472.3946 427.7087 481.7468
1976 Q3 461.7848 443.7544 479.8153 434.2097 489.3600
1976 Q4 575.3263 556.9391 593.7135 547.2055 603.4471
Resposta resumida:
AIC BIC RMSE EAM
608.7296 628.1596 12.62443 9.374296
v. Utilize os outros métodos de suavização exponencial (dentre SES, SEH e HW), medindo
também o REQM, AIC, BIC e EAM. O modelo escolhido de acordo com os testes realmente é
o mais adequado? Por quê?
Utilizando o modelo exponencial simples
> fit_ses = ses(aust, h = 4)
> summary(fit_ses)
Forecast method: Simple exponential smoothing
Model Information:
Simple exponential smoothing
AIC AICc BIC
770.6770 771.0770 777.1536
Error measures:
ME RMSE MAE MPE MAPE MASE ACF1
Training set 14.66255 49.13627 36.48202 2.502544 9.55419 2.263621 -0.06483049
3
Forecasts:
Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
1976 Q1 486.172 422.1937 550.1502 388.3257 584.0183
1976 Q2 486.172 420.5618 551.7822 385.8298 586.5141
1976 Q3 486.172 418.9694 553.3745 383.3946 588.9494
1976 Q4 486.172 417.4140 554.9300 381.0157 591.3282
Resposta resumida:
AIC BIC RMSE EAM
770.6770 777.1536 49.13627 36.48202
Utilizando o modelo de suavização exponencial de Holt.
> #suav. expo holt
> fit_seh = holt(aust, h = 4)
> summary(fit_seh)
Forecast method: Holt’s method
sigma: 43.9131
AIC AICc BIC
756.1614 757.1959 766.9558
Error measures:
ME RMSE MAE MPE MAPE MASE ACF1
Training set 4.286387 42.51872 35.42307 -0.1355001 9.584672 2.197916 -0.06359112
Forecasts:
Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
1976 Q1 503.6742 447.3972 559.9511 417.6060 589.7424
1976 Q2 509.2712 452.9791 565.5632 423.1799 595.3624
1976 Q3 514.8682 458.5422 571.1941 428.7250 601.0113
1976 Q4 520.4651 464.0789 576.8514 434.2298 606.7005
Resposta resumida:
AIC BIC RMSE EAM
756.1614 766.9558 42.51872 35.42307
vi. Com base no modelo escolhido de acordo com os testes, retorne as estimativas para os 4
trimestres de 1976. Comparando essas projeções com os valores observados reais (presentes na
série), qual seria o EQM, REQM e EAM para essas 4 observações? Opcional: faça o mesmo
para os outros dois modelos e identifique o que possui melhor generalização.
Utilizando como medidas de comparação as medidas: AIC, BIC, EQM e EAM. O modelo
escolhido foi o de Holt Winter’s aditivo.
4
aust2= ts(ausbeer[81:84],frequency=4,start=1976)
aust2
Qtr1 Qtr2 Qtr3 Qtr4
1976 510 433 453 548
valor.ajus=c(499.2172, 454.7277, 461.7848 , 575.3263)
valor.ajus
[1] 499.2172 454.7277 461.7848 575.3263
cbind(aust2,valor.ajus)
aust2 valor.ajus
1976 Q1 510 499.2172
1976 Q2 433 454.7277
1976 Q3 453 461.7848
1976 Q4 548 575.3263
> rmseadt = sqrt(sum((aust2 - valor.ajus) ^ 2) / 4)
> rmseadt
[1] 18.79003
> eamadt = ((sum(abs(aust2 - valor.ajus))) / 4)
> eamadt
[1] 17.1554
informações resumidas:
REQM EQM EAM
18.79003 353.0652 17.1554
5