Exemplo Séries Temporais: consumo Morettin e Toloi: exponential smoothing - formato ts

Licença

This work is licensed under the Creative Commons Attribution-ShareAlike 4.0 International License. To view a copy of this license, visit http://creativecommons.org/licenses/by-sa/4.0/ or send a letter to Creative Commons, PO Box 1866, Mountain View, CA 94042, USA.

License: CC BY-SA 4.0

Citação

Sugestão de citação: FIGUEIREDO, Adriano Marcos Rodrigues. Exemplo Séries Temporais: consumo Morettin e Toloi: exponential smoothing - formatos ts. Campo Grande-MS, Brasil: RStudio/Rpubs, 2023. Disponível em https://rpubs.com/amrofi/morettin_toloi_expsmooth.

1 Introdução

Os dados vem do livro de Morettin e Toloi, Análise de Séries Temporais (https://www.ime.usp.br/~pam/ST.html), consumo no varejo de São Paulo, mensais de Jan/1984 a Out/1996, em https://www.ime.usp.br/~pam/CONSUMO.XLS.

Inicialmente, trabalharemos com o formato ts, conforme o livro de Hyndman (2018), segunda edição. Depois faremos a parte da decomposição e lags com formato tsibble, conforme o livro de Hyndman (2021), terceira edição.

2 Dados em formato ts

#vou puxar os libraries que tenho costume de usar: destaque para o readxl e tseries
library(readxl) # puxar Excel
library(fpp2) # livro Hyndman formato ts - já carrega ggplot2 e forecast
library(fpp3) # livro Hyndman formato tsibble - já carrega fable, feasts e tidyr e dplyr
# codigo para puxar do Excel
# dados <- read_excel("CONSUMO morettin R.xlsx", 
#                    sheet = "dados")
dados <-structure(list(obs = structure(c(441763200, 444441600, 446947200, 
449625600, 452217600, 454896000, 457488000, 460166400, 462844800, 
465436800, 468115200, 470707200, 473385600, 476064000, 478483200, 
481161600, 483753600, 486432000, 489024000, 491702400, 494380800, 
496972800, 499651200, 502243200, 504921600, 507600000, 510019200, 
512697600, 515289600, 517968000, 520560000, 523238400, 525916800, 
528508800, 531187200, 533779200, 536457600, 539136000, 541555200, 
544233600, 546825600, 549504000, 552096000, 554774400, 557452800, 
560044800, 562723200, 565315200, 567993600, 570672000, 573177600, 
575856000, 578448000, 581126400, 583718400, 586396800, 589075200, 
591667200, 594345600, 596937600, 599616000, 602294400, 604713600, 
607392000, 609984000, 612662400, 615254400, 617932800, 620611200, 
623203200, 625881600, 628473600, 631152000, 633830400, 636249600, 
638928000, 641520000, 644198400, 646790400, 649468800, 652147200, 
654739200, 657417600, 660009600, 662688000, 665366400, 667785600, 
670464000, 673056000, 675734400, 678326400, 681004800, 683683200, 
686275200, 688953600, 691545600, 694224000, 696902400, 699408000, 
702086400, 704678400, 707356800, 709948800, 712627200, 715305600, 
717897600, 720576000, 723168000, 725846400, 728524800, 730944000, 
733622400, 736214400, 738892800, 741484800, 744163200, 746841600, 
749433600, 752112000, 754704000, 757382400, 760060800, 762480000, 
765158400, 767750400, 770428800, 773020800, 775699200, 778377600, 
780969600, 783648000, 786240000, 788918400, 791596800, 794016000, 
796694400, 799286400, 801964800, 804556800, 807235200, 809913600, 
812505600, 815184000, 817776000, 820454400, 823132800, 825638400, 
828316800, 830908800, 833587200, 836179200, 838857600, 841536000, 
844128000), class = c("POSIXct", "POSIXt"), tzone = "UTC"), t = c(1, 
2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 
20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 
36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 
52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 
68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 
84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 
100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 
113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 
126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 
139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 
152, 153, 154), consumo = c(114.13, 110.79, 116.46, 111.57, 120.66, 
121.15, 121.27, 127.02, 129.04, 133.3, 130.6, 179.39, 120.64, 
114.05, 130.6, 118.26, 145.54, 135.13, 153.35, 159.95, 150.01, 
164.93, 170.37, 220.96, 134.26, 133.11, 147.84, 164.46, 181.86, 
170.44, 186.64, 174.21, 181.62, 194.16, 181.9, 232.01, 140.16, 
130.78, 119.04, 120.73, 129.81, 111.04, 122.75, 133.95, 125.41, 
132.05, 129.54, 176.37, 110.09, 113.25, 124.03, 110.63, 116.72, 
124.63, 124.38, 130.27, 119.87, 115.75, 122.44, 162.43, 105.89, 
115.59, 147, 131.7, 131.32, 136.66, 126.43, 134.88, 128.26, 125.32, 
124.61, 166.11, 116.25, 96.93, 89.27, 101.87, 125.57, 113.31, 
109.39, 127.33, 120.56, 117.73, 113.81, 147.25, 100.15, 95.11, 
112.26, 109.39, 114.2, 113.8, 126.47, 128.36, 115.71, 116.09, 
99.53, 127.27, 87.08, 85.67, 82.02, 98.2, 96.44, 90.23, 97.15, 
95.08, 94, 93, 96.09, 129.21, 75.39, 77.7, 97.34, 84.97, 87.55, 
86.64, 90.52, 95.4, 95.2, 95.8, 101.23, 128.49, 85.63, 82.77, 
96.55, 81.33, 96.91, 83.76, 90.19, 114.84, 108.4, 106.05, 109.71, 
143.86, 99.12, 99.28, 114.75, 106.13, 110.02, 108.07, 112.52, 
113.87, 107.84, 112.12, 112.03, 139.37, 92.24, 93.56, 107.37, 
102.89, 114.78, 102.88, 118.41, 119.23, 117.36, 122.06)), row.names = c(NA, 
-154L), class = c("tbl_df", "tbl", "data.frame"))
# atribuir objeto ts
dados.ts<- ts(dados$consumo,start=c(1984,1), frequency = 12)
# plot inicial
plot(dados.ts,main="Consumo do varejo de Sao Paulo, Morettin e Toloi (2006)",
              xlab="Ano",ylab="Indice")

#
# plot no ggplot
library(ggplot2)
datas <- seq(as.Date(paste(c(start(dados.ts),1), collapse="/")), 
             by = "month", length.out = length(dados$consumo))
dados.df <- data.frame(date = datas, value = dados$consumo)
ggplot(dados.df) +
  aes(x = date, y = value) +
  geom_line(colour = "#112446") +
  labs(
    x = "Mês/Ano",
    y = "Índice (base 100)",
    title = "Série de consumo:",
    subtitle = "varejo de São Paulo",
    caption = "Fonte: Morettin e Toloi (2006). Disponível em: <https://www.ime.usp.br/~pam/CONSUMO.XLS>"
  ) +
  ggthemes::theme_economist_white()

2.0.1 Estatísticas descritivas

Olhando as estatísticas descritivas, o pesquisador percebe a amplitude de variação, e o gráfico do dygraphs permite “caminhar” com o mouse sobre a série. Perceberás que a sazonalidade está bem marcada com picos em dezembro de cada ano.

# estatisticas basicas
summary(dados.ts)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   75.39  102.12  116.36  120.94  130.52  232.01

#   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
#  75.39  102.12  116.36  120.94  130.52  232.01 
# pelo pacote dygraph dá mais opções
library(dygraphs)
dygraph(dados.ts)

É possível visualizar nos plots acima: sazonalidade (por exemplo, picos em dezembro de cada ano); a tendência aparentemente crescente e depois decresce com a “crise” brasileira; e uma aparente não-estacionariedade (média e variância mudam no tempo). Em outra postagem aplicarei o teste de raiz unitária na série para avaliar a estacionariedade de modo mais explícito.

3 Suavização (ou alisamento) Exponencial Simples (SES) ou ainda, ‘Holt 1 parâmetro’

Este método é apropriado para séries que se movem aleatoriamente em torno de uma média ou nível constante, sem tendência e sem sazonalidade: \(Z_t=μ_t+a_t\) em que \(a_t\) é um ruído branco de média zero e variância constante e \(μ_t\) é um parâmetro localmente constante (nível); \(Z_t\) é a série em análise.

\[ {\bar Z_t} = \alpha {Z_t} + \left( {1 - \alpha } \right){\bar Z_{t - 1}} \]

Em que α é a constante de suavização, 0<α<1. Portanto, similarmente,

\[ {\bar Z_{t - 1}} = \alpha {Z_{t - 1}} + \left( {1 - \alpha } \right){\bar Z_{t - 2}} \]

O que substituindo (2) em (1) e desenvolvendo os termos tem-se:

\[ {\bar Z_t} = \alpha {Z_t} + \alpha \left( {1 - \alpha } \right){Z_{t - 1}} + \alpha {\left( {1 - \alpha } \right)^2}{Z_{t - 2}} + \ldots \]

Ou seja, a SES faz uma previsão que é a média ponderada por \(α(1-α)^k\) que, quanto mais recente, maior é o peso.

Se \(0<α<1\), então \(α>α(1-α)>α(1-α)^2>\ldots>α(1-α)^k\)

Na MMS tem-se o mesmo peso entre observações recentes e antigas. Na SES, a vantagem sobre MMS é que aplica pesos maiores às observações recentes. A questão será achar α que minimiza a soma de quadrados de ajustamento:

\[ S = \mathop \sum \limits_{t = l + 1}^N {\left( {{Z_t} - {{\hat Z}_{t - 1}}} \right)^2} \]

O forecast será constante para todas as observações futuras!!!

3.1 Código SES

x<-dados.ts  # coringa para automatizar
# simple exponential smoothing (SES)

x.ses<-ses(x,h=12,level = c(95))
x.ses          # exibe o forecast por holt 1 parametro

##          Point Forecast    Lo 95    Hi 95
## Nov 1996       117.0958 80.15345 154.0382
## Dec 1996       117.0958 78.19363 155.9980
## Jan 1997       117.0958 76.32792 157.8637
## Feb 1997       117.0958 74.54393 159.6477
## Mar 1997       117.0958 72.83179 161.3598
## Apr 1997       117.0958 71.18345 163.0082
## May 1997       117.0958 69.59227 164.5994
## Jun 1997       117.0958 68.05269 166.1389
## Jul 1997       117.0958 66.55999 167.6316
## Aug 1997       117.0958 65.11014 169.0815
## Sep 1997       117.0958 63.69964 170.4920
## Oct 1997       117.0958 62.32545 171.8662

summary(x.ses)  # fornece o relatorio da estimacao

## 
## Forecast method: Simple exponential smoothing
## 
## Model Information:
## Simple exponential smoothing 
## 
## Call:
##  ses(y = x, h = 12, level = c(95)) 
## 
##   Smoothing parameters:
##     alpha = 0.33 
## 
##   Initial states:
##     l = 115.5709 
## 
##   sigma:  18.8485
## 
##      AIC     AICc      BIC 
## 1684.099 1684.259 1693.210 
## 
## Error measures:
##                      ME    RMSE      MAE       MPE     MAPE     MASE       ACF1
## Training set 0.03000474 18.7257 13.03585 -1.699482 10.66839 0.813645 0.05177914
## 
## Forecasts:
##          Point Forecast    Lo 95    Hi 95
## Nov 1996       117.0958 80.15345 154.0382
## Dec 1996       117.0958 78.19363 155.9980
## Jan 1997       117.0958 76.32792 157.8637
## Feb 1997       117.0958 74.54393 159.6477
## Mar 1997       117.0958 72.83179 161.3598
## Apr 1997       117.0958 71.18345 163.0082
## May 1997       117.0958 69.59227 164.5994
## Jun 1997       117.0958 68.05269 166.1389
## Jul 1997       117.0958 66.55999 167.6316
## Aug 1997       117.0958 65.11014 169.0815
## Sep 1997       117.0958 63.69964 170.4920
## Oct 1997       117.0958 62.32545 171.8662

plot(x.ses,col=1, 
     main = "Índice de volume de vendas no varejo Total 
     de São Paulo", sub="modelo SES")

Atentar que é linear porque não tem tendência, observada no summary.

3.2 Código Holt

# Holt dois parametros

x.holt<-holt(x,h=12,level = c(95))
x.holt          # exibe o forecast por holt 1 parametro

##          Point Forecast    Lo 95    Hi 95
## Nov 1996       117.2237 80.01095 154.4364
## Dec 1996       117.2512 78.02936 156.4731
## Jan 1997       117.2788 76.14468 158.4129
## Feb 1997       117.3064 74.34396 160.2688
## Mar 1997       117.3339 72.61689 162.0509
## Apr 1997       117.3615 70.95512 163.7678
## May 1997       117.3890 69.35177 165.4263
## Jun 1997       117.4166 67.80106 167.0321
## Jul 1997       117.4442 66.29814 168.5902
## Aug 1997       117.4717 64.83882 170.1046
## Sep 1997       117.4993 63.41951 171.5790
## Oct 1997       117.5268 62.03708 173.0166

summary(x.holt)  # fornece o relatorio da estimacao

## 
## Forecast method: Holt's method
## 
## Model Information:
## Holt's method 
## 
## Call:
##  holt(y = x, h = 12, level = c(95)) 
## 
##   Smoothing parameters:
##     alpha = 0.3329 
##     beta  = 1e-04 
## 
##   Initial states:
##     l = 109.859 
##     b = 0.0267 
## 
##   sigma:  18.9864
## 
##      AIC     AICc      BIC 
## 1688.305 1688.711 1703.490 
## 
## Error measures:
##                      ME     RMSE      MAE       MPE     MAPE      MASE
## Training set 0.05516869 18.73824 13.04413 -1.672203 10.68162 0.8141618
##                    ACF1
## Training set 0.05016778
## 
## Forecasts:
##          Point Forecast    Lo 95    Hi 95
## Nov 1996       117.2237 80.01095 154.4364
## Dec 1996       117.2512 78.02936 156.4731
## Jan 1997       117.2788 76.14468 158.4129
## Feb 1997       117.3064 74.34396 160.2688
## Mar 1997       117.3339 72.61689 162.0509
## Apr 1997       117.3615 70.95512 163.7678
## May 1997       117.3890 69.35177 165.4263
## Jun 1997       117.4166 67.80106 167.0321
## Jul 1997       117.4442 66.29814 168.5902
## Aug 1997       117.4717 64.83882 170.1046
## Sep 1997       117.4993 63.41951 171.5790
## Oct 1997       117.5268 62.03708 173.0166

plot(x.holt,col=1, 
     main = "Índice de volume de vendas no varejo Total 
     de São Paulo", sub="modelo Holt")

Atentar que parece linear mas tem uma pequena tendência que pode ser observada no summary.

3.3 Código Holt-Winters - 3 parâmetros

# Holt-Winters 
x.hwm <- hw(x, seasonal="multiplicative")  # multiplicativo
plot(x,ylab="Indice base 100", main=
       "Índice de volume de vendas no varejo Total de São Paulo",
     type="o",  xlab="Mês/Ano")
lines(fitted(x.ses), col="red", lty=2)
lines(x.ses$mean, type="o", col="red")
lines(x.hwm$mean, type="o", col="green")
legend("topright",lty=1, pch=1, col=1:3, 
       c("original","Holt","HW Multiplicativo"))

summary(x.hwm)

## 
## Forecast method: Holt-Winters' multiplicative method
## 
## Model Information:
## Holt-Winters' multiplicative method 
## 
## Call:
##  hw(y = x, seasonal = "multiplicative") 
## 
##   Smoothing parameters:
##     alpha = 0.5328 
##     beta  = 0.0102 
##     gamma = 1e-04 
## 
##   Initial states:
##     l = 125.3108 
##     b = 1.595 
##     s = 1.3325 1.0103 1.0256 1.0065 1.0415 1.0014
##            0.9618 1.0115 0.9279 0.9403 0.8566 0.884
## 
##   sigma:  0.0707
## 
##      AIC     AICc      BIC 
## 1448.858 1453.358 1500.486 
## 
## Error measures:
##                      ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -0.8609658 8.202305 5.964579 -0.9323561 5.014909 0.3722848
##                   ACF1
## Training set 0.2406417
## 
## Forecasts:
##          Point Forecast     Lo 80    Hi 80     Lo 95    Hi 95
## Nov 1996       118.9314 108.15746 129.7054 102.45407 135.4088
## Dec 1996       157.1159 140.91698 173.3148 132.34179 181.8900
## Jan 1997       104.4061  92.43792 116.3744  86.10232 122.7100
## Feb 1997       101.3382  88.62472 114.0516  81.89463 120.7817
## Mar 1997       111.4200  96.29497 126.5451  88.28824 134.5518
## Apr 1997       110.1296  94.09179 126.1674  85.60189 134.6573
## May 1997       120.2463 101.58710 138.9054  91.70953 148.7830
## Jun 1997       114.5304  95.69544 133.3654  85.72481 143.3360
## Jul 1997       119.4445  98.71949 140.1694  87.74835 151.1406
## Aug 1997       124.4382 101.74268 147.1338  89.72838 159.1481
## Sep 1997       120.4505  97.43228 143.4688  85.24716 155.6539
## Oct 1997       122.9276  98.38061 147.4746  85.38621 160.4690
## Nov 1997       121.3019  96.05134 146.5524  82.68452 159.9192
## Dec 1997       160.2422 125.54335 194.9411 107.17488 213.3095
## Jan 1998       106.4802  82.53968 130.4207  69.86633 143.0941
## Feb 1998       103.3479  79.26163 127.4342  66.51111 140.1848
## Mar 1998       113.6261  86.21665 141.0356  71.70696 155.5453
## Apr 1998       112.3065  84.30425 140.3088  69.48074 155.1323
## May 1998       122.6193  91.05624 154.1823  74.34778 170.8908
## Jun 1998       116.7869  85.78721 147.7866  69.37697 164.1968
## Jul 1998       121.7939  88.49089 155.0970  70.86133 172.7265
## Aug 1998       126.8819  91.17601 162.5878  72.27444 181.4894
## Sep 1998       122.8120  87.27497 158.3491  68.46279 177.1613
## Oct 1998       125.3338  88.07275 162.5948  68.34796 182.3195

É possível perceber que o Holt-Winters (HW - 3 parâmetros) capta melhor a sazonalidade, que é justamente a diferença deste comparado ao Holt (2 parâmetros).

4 Avaliação da simulação

Para uma só equação, podem-se olhar os indicadores tradicionais como o R², estatísticas F e t, critérios de informação de Akaike, Schwarz entre outros. Para modelos de previsão, devo ter um erro de previsão mínimo, definido pelo pesquisador de modo a julgá-lo adequado.

As vezes pode-se ter equações bem estimadas junto de outras com poucos parâmetros significativos, num mesmo sistema. Dependendo do objetivo a ser alcançado, se simular ou explicar, busca-se ajustamento com alta significância (para explicar) ou busca-se melhor poder de explicação. Por melhor poder de explicação entende-se: valores previstos mais próximos dos observados.

Alguns critérios para avaliação da previsão podem ser citados: Raiz do Erro de Simulação Quadrático Médio; Raiz do Erro Quadrático Médio Percentual; Erro de Simulação Médio; Erro Percentual Médio; Erro Absoluto Médio; Erro Absoluto Percentual Médio; Raiz do Erro de Simulação Ex-post Quadrático Médio; Coeficiente de desigualdade de Theil (U).

4.1 Raiz do Erro de Simulação Quadrático Médio (RMSE)

Este critério é dado pela expressão:

\[ RMSE = \;\sqrt {\frac{1}{T}\mathop \sum \limits_{t = 1}^T {{\left( {Y_t^s - Y_t^a} \right)}^2}} \]

Em que: \(Y_t^s\) é o valor simulado (previsto) de \(Y_t\); \(Y_t^a\) é o valor efetivo (observado); e T é o número de períodos na simulação.

O resultado é um escalar com o valor da raiz quadrada da média da diferença quadrática entre x e y, ou seja, a raiz do erro quadrático médio. Normalmente deve ser calculado pela fórmula, mas é reportado pelo programa quando se faz um “forecast”.

Para o cálculo da Raiz do Erro de Simulação Ex-post Quadrático Médio, utilizam-se os resultados do modelo de previsão aplicado a um conjunto restrito (histórico) para checar os erros no período de simulação. Exemplo: usando os dados de 1950:01 a 1985:04, históricos, simulam-se os resultados para 1986:01 a 1988:01 e calcula-se o erro ex-post.

4.2 Raiz do Erro Quadrático Médio Percentual (PRMSE)

\[ PRMSE = \;\sqrt {\frac{1}{T}\mathop \sum \limits_{t = 1}^T {{\left( {\frac{{Y_t^s - Y_t^a}}{{Y_t^a}}} \right)}^2}} \]

Em que: \(Y_t^s\) é o valor simulado (previsto) de \(Y_t\); \(Y_t^a\) é o valor efetivo (observado); e T é o número de períodos na simulação. Neste caso, a diferença entre o previsto e o observado (o erro) é dividida pelo valor observado, obtendo-se assim em base percentual relativo ao valor observado.

4.3 Erro de Simulação Médio (ME) e o Erro Percentual Médio (MPE)

O erro de simulação é basicamente o termo \((Y_t^s - Y_t^a)\), como definido anteriormente. O erro médio (ME) é, portanto, a média dos erros para T observações, enquanto o erro percentual médio (MPE) é obtido relativizando-se pelos valores observados, como nas expressões a seguir:

\[ ME= \frac{1}{T}\sum\limits_{t = 1}^T {\left( {Y_t^s - Y_t^a} \right)} \]

\[ MPE = 100.\frac{1}{T}\sum\limits_{t = 1}^T {\frac{{\left( {Y_t^s - Y_t^a} \right)}}{{Y_t^a}}} \]

Uma ressalva necessária é que: “[…] O problema com os erros médios é que eles podem ser quase zero, quando elevados erros positivos cancelam elevados erros negativos” (PINDYCK e RUBINFELD, 2004, p.447).

Alguns autores utilizam estas duas últimas medidas com o termo entre parênteses $ {( {Y_t^s - Y_t^a} )} $ na forma de ‘valor absoluto’ \(\left| {\left( {Y_t^s - Y_t^a} \right)} \right|\) . Tal procedimento auxilia a investigar o erro de previsão do modelo evitando os problemas de ‘erros positivos’ cancelarem ‘erros negativos’. Desta forma, surgem as outras duas opções: MAE e MAPE.

4.4 Erro Absoluto Médio (MAE) e Erro Absoluto Percentual Médio (MAPE)

O erro absoluto médio (MAE) considera a expressão análoga do erro de simulação médio com os valores absolutos como na expressão:

\[ MAE = \frac{1}{T}\sum\limits_{t = 1}^T {\left| {Y_t^s - Y_t^a} \right|} \]

Similarmente ao erro percentual médio, o erro absoluto percentual médio (MAPE) será:

\[ MAPE = 100.\frac{1}{T}\sum\limits_{t = 1}^T {\left| {\frac{{\left( {Y_t^s - Y_t^a} \right)}}{{Y_t^a}}} \right|} \]

Desejam-se valores baixos das medidas. Entretanto, é importante observar se o modelo está bem ajustado para os ‘pontos de inflexão’ nos dados (oscilações, vales e picos que se observam nos gráficos previstos e realizados).

4.5 Erro Absoluto Percentual Médio Simétrico (sMAPE)

As medidas de acurácia tipo MAPE têm a desvantagem de penalizar mais os erros negativos que os positivos. Assim, surge a medida sMAPE ou MAPE “simétrico” proposto por Armstrong (1978 , p. 348), citado por Hyndman e Athanasopoulos (2018). Esta medida foi utilizada para escolha do melhor modelo na competição de previsão do M3 de Makridakis (\(M^3\) competition). A expressão para cálculo é do tipo:

\[ \text{sMAPE} = \text{mean}\left(200|{Y_t^s - Y_t^a}|/(Y_t^s+Y_t^a)\right) \]

Esta medida não está disponibilizada diretamente nas saídas do pacote forecast nem no fpp2, em virtude da crítica de Hyndman e Athanasopoulos (2018) em que citam Hyndman e Koehler (2006). A medida sMAPE pode ser obtida pelos pacotes Metrics ou TSPred.

4.6 Erro absoluto médio escalonado (MASE)

Para séries temporais sazonais, um erro escalonado para a sazonalidade pode ser definido usando previsões ingênuas sazonais (Hyndman e Athanasopoulos, 2021):

\[ q_{j} = \frac{\displaystyle e_{j}} {\displaystyle\frac{1}{T-m}\sum_{t=m+1}^T |y_{t}-y_{t-m}|}. \]

A métrica MASE será então dada pela expressão:

$$ = (|q_{j}|).

$$

4.7 critérios de informação, Akaike (AIC), Schwarz (SIC ou BIC) e Hannan-Quin (HQ)

A ideia é que nos critérios de informação, acrescentando lags se reduz a soma dos resíduos (pois mais variação é explicada no modelo), ou seja, foi acrescentada mais informação. Mas ainda deve-se penalizar a adição do regressor, uma ideia de redução dos graus de liberdade do modelo.

Em geral, os critérios têm a seguinte expressão: \(C=lnσ ̂^2 (T)+c_T φ(T)\).

O primeiro termo associa a variância, e seu valor deve reduzir com o aumento de variáveis. No segundo termo tem-se: \(c_T\) é o número de parâmetros estimados; e \(φ(T)\) é a ordem do processo que penaliza o acréscimo de variáveis. Assim, quanto menor C melhor é o modelo, ou seja, tem menos resíduos, menos parte não explicada. O Critério de Informação de Akaike (ou AIC de Akaike’s Information Criterion) ou o Critério de Informação de Schwarz ou Bayesiano (ou SIC de Schwarz’s Information Criterion ou em alguns livros BIC de Bayesian Information Criterion) são expressos no Eviews da forma já logaritmizada como:

\[AIC = - \frac{{2l}}{n} + \frac{{2k}}{n}\]

\[SIC = - \frac{{2l}}{n} + \frac{{k.\log n}}{n}\]

\[ l = - {\frac{n}{2}} {\left( {1 + \log(2\pi)+\log \left( {\frac{{\hat \varepsilon ^\prime \hat \varepsilon }}{n}} \right)} \right)} \]

Em que: k é o número de regressores incluindo-se o intercepto; n é o número de observações; l é o log Verossimilhança da regressão; e \(\hat \varepsilon\) são os resíduos estimados do modelo. No formato mais simplificado exposto por Greene (2002), tem-se igual ao de Bueno (2008, p.47):

\[ \begin{array}{l} SIC = BIC\left( {p,q} \right) = \ln {{\hat \sigma }^2} + k\frac{{\ln T}}{T}\\ {{\hat \sigma }^2} = \frac{{\sum\nolimits_{t = 1}^T {\hat \varepsilon _t^2} }}{T} \end{array} \]

Em que k=p+q+1 se houver intercepto.

\[ AIC\left( {p,q} \right) = \ln {{\hat \sigma }^2} + k\frac{{2}}{T} \]

\[ HQ\left( {p,q} \right) = \ln {{\hat \sigma }^2} + k\frac{2}{T} \ln T \]

Quanto mais parâmetros são estimados, menor será o erro estimado, mas isto é penalizado no segundo termo dos critérios. Desejo menores valores dos critérios. O BIC é consistente assintoticamente, e o AIC funciona melhor em pequenas amostras, mas é pior com mais parâmetros. Já o HQ é menos forte que o BIC, ou seja, o preferido é o BIC.

Podemos resumir os resultados de acurácia para visualizar que o HW realmente foi o melhor modelo em termos de AICc, o menor, e nas métricas de acurácia, também o menor principalmente olhando o RMSE e MASE (Figura 1).

Figura 1. Resumo dos resultados dos critérios de informação e das métricas de acurácia, modelos SES, Holt e Holt-Winters.

Fonte: Elaboração própria.

5 Referências

HYNDMAN, Rob. (2018). fpp2: Data for “Forecasting: Principles and Practice” (2nd Edition). R package version 2.3. Disponível em: https://CRAN.R-project.org/package=fpp2.

HYNDMAN, R.J., & ATHANASOPOULOS, G. (2018) Forecasting: principles and practice, 2nd edition, OTexts: Melbourne, Australia. Disponível em: https://otexts.com/fpp2/. Accessed on 27 Mar 2023.

HYNDMAN, R.J., & ATHANASOPOULOS, G. (2021) Forecasting: principles and practice, 3rd edition, OTexts: Melbourne, Australia. Disponível em: https://otexts.com/fpp3/.

MORETTIN, Pedro A.; TOLOI, Clélia M.C. Análise de Séries Temporais. São paulo: Editora Edgard Blücher, 2004. https://www.ime.usp.br/~pam/st/