Tag

tsibble Archives - Análise Macro

{tidyverts}: séries temporais no R

By | Data Science

{tidyverts} é uma família de pacotes de R criada para ser a próxima geração de ferramentas para modelagem e previsão de séries temporais, substituindo o famoso pacote {forecast}.

Utilizando uma interface simples e integrada com os pacotes do tidyverse, é possível construir uma ampla gama de modelos de previsão univariados e multivariados: ARIMA, VAR, suavização exponencial via espaço de estado (ETS), modelo linear (TSLM), autorregressivo (AR), passeio aleatório (RW), autoregressão de rede neural (NNETAR), Prophet, etc.

Neste exercício, daremos uma breve introdução aplicada a modelagem e previsão de séries temporais usando os pacotes do {tidyverts}.

Conhecendo os pacotes

Um breve resumo do que a família de pacotes do {tidyverts} tem a oferecer:

fable

  • Coleção de modelos univariados e multivariados de previsão
  • Modelagem de séries temporais em formato "tidy"
  • Especificação de modelos utiliza terminologia de fórmula (y ~ x)

fabletools

  • Extensões e ferramentas para construção de modelos
  • Combinação de modelos, previsão hierárquica e extração de resultados
  • Obtenção de medidas de acurácia e visualização de dados

feasts

  • Decomposição de séries temporais
  • Extração e visualização de componentes de séries temporais
  • Análise de autocorrelação, testes de raiz unitária, etc.

tsibble

  • Estrutura de dados tidy para séries temporais no R
  • Funções para tratamento de dados
  • Objeto orientado aos dados e a modelos, integrado ao tidyverse

Fluxo de trabalho para previsão

Com o tidyverts o processo de construir um modelo de previsão pode ser dividido em poucos passos:

Com esse esquema em mente, vamos ilustrar esse processo com um exercício prático e didático: construir um modelo de previsão para a taxa de crescimento do PIB brasileiro.

Pacotes

Para reproduzir o exercício a seguir você precisará dos seguintes pacotes:

Dados tidy

Utilizaremos o dataset global_economy armazenado como um objeto tsibble, trazendo variáveis econômicas em frequência anual para diversos países. Nosso interesse é a série da taxa de crescimento do PIB brasileiro:

Visualização de dados

Visualização é uma etapa essencial para entender os dados, o que permite identificar padrões e modelos apropriados. No nosso exemplo, criamos um gráfico de linha para plotar a série do PIB brasileiro usando a função autoplot():

Podemos também plotar os correlogramas ACF e PACF para identificar o processo estocástico da série, obtendo alguns modelos candidatos:

Especificação do modelo

Existem muitos modelos de séries temporais diferentes que podem ser usados para previsão, e especificar um modelo apropriado para os dados é essencial para produzir previsões.

Os modelos no framework do fable são especificados usando funções com nomenclatura abreviada do nome do modelo (por exemplo, ARIMA(), AR(), VAR(), etc.), cada uma usando uma interface de fórmula (y ~ x). As variáveis de resposta são especificadas à esquerda da fórmula e a estrutura do modelo é escrita à direita.

Por exemplo, um modelo ARIMA(1,0,2) para a taxa de crescimento do PIB pode ser especificado com: ARIMA(Growth ~ pdq(1, 0, 2)).

Neste caso, a variável resposta é Growth e está sendo modelada usando a estrutura de um modelo ARMA(1, 2) especificada na função especial pdq().

Existem diversas funções especiais para definir a estrutura do modelo e em ambos os lados da fórmula pode ser aplicado transformações. Consulte detalhes da documentação do fable.

Estimar o modelo

Identificado um modelo (ou mais) apropriado, podemos em seguida fazer a estimação usando a função model()1.Neste exemplo, estimaremos os seguintes modelos: ARIMA(1,0,2), ARIMA(1,0,0), ARIMA(0,0,2), o algoritmo de seleção automatizada do auto ARIMA criado pelo prof. Rob Hyndman e um passeio aleatório.

Diagnóstico do modelo

O objeto resultante é uma "tabela de modelo" ou mable, com a saída de cada modelo em cada coluna:

Para obter os critérios de informação use a função glance():

Os critérios de informação indicam que, dos modelos estimados, o modelo automatizado ARIMA(1,1,1) apresentou o menor valor de AICc - seguido pelos demais identificados pelos correlogramas ACF e PACF. Com a função gg_tsresiduals() podemos verificar o comportamento dos resíduos deste modelo, indicando que os resíduos se comportam como ruído branco:

Um teste de autocorrelação (Ljung Box) retorna um p-valor grande, também indicando que os resíduos são ruído branco:

Também pode ser interessante visualizar o ajuste do modelo. Utilize a função augment() para obter os valores estimados:

Previsão

Com o modelo escolhido, previsões podem ser geradas com a função forecast() indicando um horizonte de escolha.

Perceba que os pontos de previsão médios gerados são bastante similares a um processo de passeio aleatório (equivalente a um ARIMA(0,1,0)). O trabalho adicional de especificar termos AR e MA trouxe pouca diferença para os pontos de previsão neste exemplo, apesar de ser perceptível que os intervalos de confiança do modelo auto ARIMA são mais estreitos do que de um passeio aleatório.

Além disso, a previsão fora da amostra gerada ficou bastante aquém dos dados reais para a taxa de crescimento do PIB brasileiro observados no horizonte em questão, configurando apenas um exercício didático.

Saiba mais

Estes são apenas alguns dos recursos e ferramentas disponíveis na família de pacotes do tidyverts. Para uma referência aprofundada, confira o livro Forecasting: Principles and Practice, 3rd Edition, de Hyndman e Athanasopoulos (2021).

Aprenda com maior profundidade sobre sobre os temas abordados nos cursos de Séries Temporais e de Previsão Macro. Faça parte do Clube AM para acesso completo aos códigos de R e Python deste e de outros exercícios.

Confira outros exercícios aplicados com pacotes do tidyverts:

 


[1] A função suporta estimação dos modelos com computação paralela usando o pacote future, veja detalhes na documentação e este post para saber mais sobre o tema.

 

 

Fazendo limpeza de dados Macroeconômicos com R e o Python

By | Hackeando o R, Python

O processo de limpeza de dados é extremamente importante para a análise de dados. Sem esta etapa, seria impossível lidar com dados "sujos", dos quais possuem formatos, classes e estruturas não compatíveis com o R e o Python, principalmente para o modo de lidar com dados Macroeconômicos. O objetivo do post de hoje será mostrar os pontos principais que devem ser levados em consideração ao realizar o processo de limpeza de dados.

Ao realizar a importação de dados Macroeconômicos, é comum enfrentar alguns desafios para efetivamente lidar com o dataset, sendo eles:

  • Nomes de colunas indesejadas;
  • Colunas (variáveis) não necessárias;
  • Filtragem;
  • Formatos/tipos dos dados de variáveis não corretos;
  • Problemas com o formato de data

Cada problema listado acima pode ser resolvido facilmente utilizando o R e o Python. Elencaremos cada problema em ambas as linguagens

Limpeza de dados Macroeconômicos no R

Tomaremos como exemplo a série do PIB a preços de mercado em milhões retirado do Ipeadata utilizando o pacote {ipeadatar}. O primeiro passo será investigar o data frame e sua estrutura

Vemos que o data frame já é importado em classe tibble com cinco colunas, representando o código, a data, o valor da série, o nome do território e o código territorial, respectivamente. Não há motivos para manter as colunas que não sejam a de data e do valor da série, portanto, iremos remover estas colunas utilizando a função select().

O próximo passo, este mais a gosto do usuário, será renomear as colunas do objeto. Como queremos que sejam facilmente reconhecidos as variáveis, utilizaremos nomes que sugerem o seu significado.

Um ponto interessante do R é a facilidade de manipulação com séries temporais. Ao inspecionar o tibble, vemos que há um intervalo trimestral para as observações, seria interessante, portanto, utilizar um formato de data mais conveniente para o caso. Podemos realizar a mudança utilizando a função yearquarter() do pacote {tsibble}.

Fácil até então, certo? Porém, nem sempre pode ser tão facilitado assim. Um caso em que as variáveis não são importadas com seus devidos formatos ocorre ao utilizar o pacote {sidrar} para obter dados do Sidra - IBGE.

Como exemplo podemos utilizar a Pesquisa Mensal de Serviços da Tabela 6442. Obtemos sua API através do site do Sidra.

Ao analisar o data frame importado, não parece tão simples, certo? O primeiro passo para a limpeza, será definir quais variáveis necessitamos e destas, quais devem ser realizadas as manipulações. Para o momento, queremos o valor da variável, a coluna que define a categoria que o valor se encaixa, uma representado se o valor possui ou não ajuste sazonal, uma coluna definindo qual medida se refere a variável, sendo os valores a receita nominal ou o volume, e não menos importante, a data das observações.

O processo será selecionar estas colunas citadas, mudar seus nomes, filtrar a categoria de interesse da medida, mudar de long para wide e por fim, realizar as mudanças necessárias para a data.

Limpeza de dados Macroeconômicos no Python

Todos os procedimentos acima podem ser feitos também utilizando o Python. Seremos mais diretos desta vez, importando a tabela 1442 do Sidra com a biblioteca {sidrapy} que se refere aos dados do IPCA.

O primeiro procedimento será avaliar  o data frame importado e realizar os ajustes necessários.

Veja que os dados importados possuem colunas com nomes diferentes daqueles importados pelo R, a diferença estão que as colunas com os verdadeiros nomes estão localizados na primeiro observação do data frame, ou seja, precisamos substituir as colunas pela primeira observação.

O segundo passo será selecionar apenas as variáveis de interesse, bem como alterar os seus nomes para algo mais útil.

Por fim, devemos realizar o filtro para a variação do IPCA mensal, de acordo com a categoria inserida na coluna 'variable' e em seguida, alterar a coluna de data para o tipo datetime e a realocar para o índice do data frame.

Quer saber mais sobre R e Python?

Veja nossos cursos de:

{tidyverts}: séries temporais no R

By | Data Science

{tidyverts} é uma família de pacotes de R criada para ser a próxima geração de ferramentas para modelagem e previsão de séries temporais, substituindo o famoso pacote {forecast}.

Utilizando uma interface simples e integrada com os pacotes do tidyverse, é possível construir uma ampla gama de modelos de previsão univariados e multivariados: ARIMA, VAR, suavização exponencial via espaço de estado (ETS), modelo linear (TSLM), autorregressivo (AR), passeio aleatório (RW), autoregressão de rede neural (NNETAR), Prophet, etc.

Neste exercício, daremos uma breve introdução aplicada a modelagem e previsão de séries temporais usando os pacotes do {tidyverts}.

Conhecendo os pacotes

Um breve resumo do que a família de pacotes do {tidyverts} tem a oferecer:

fable

  • Coleção de modelos univariados e multivariados de previsão
  • Modelagem de séries temporais em formato "tidy"
  • Especificação de modelos utiliza terminologia de fórmula (y ~ x)

fabletools

  • Extensões e ferramentas para construção de modelos
  • Combinação de modelos, previsão hierárquica e extração de resultados
  • Obtenção de medidas de acurácia e visualização de dados

feasts

  • Decomposição de séries temporais
  • Extração e visualização de componentes de séries temporais
  • Análise de autocorrelação, testes de raiz unitária, etc.

tsibble

  • Estrutura de dados tidy para séries temporais no R
  • Funções para tratamento de dados
  • Objeto orientado aos dados e a modelos, integrado ao tidyverse

Fluxo de trabalho para previsão

Com o tidyverts o processo de construir um modelo de previsão pode ser dividido em poucos passos:

Com esse esquema em mente, vamos ilustrar esse processo com um exercício prático e didático: construir um modelo de previsão para a taxa de crescimento do PIB brasileiro.

Pacotes

Para reproduzir o exercício a seguir você precisará dos seguintes pacotes:

library(fable) # CRAN v0.3.1
library(fabletools) # CRAN v0.3.2
library(tsibble) # CRAN v1.1.1
library(tsibbledata) # CRAN v0.4.0
library(feasts) # CRAN v0.2.2
library(dplyr) # CRAN v1.0.7
library(tidyr) # CRAN v1.2.0
library(ggplot2) # CRAN v3.3.5

 

Dados tidy

Utilizaremos o dataset global_economy armazenado como um objeto tsibble, trazendo variáveis econômicas em frequência anual para diversos países. Nosso interesse é a série da taxa de crescimento do PIB brasileiro:

pib_br <- tsibbledata::global_economy %>% 
dplyr::filter(Country == "Brazil") %>% 
dplyr::select(Year, Growth) %>% 
tidyr::drop_na()
pib_br

# # A tsibble: 57 x 2 [1Y]
# Year Growth
# <dbl> <dbl>
# 1 1961 10.3 
# 2 1962 5.22 
# 3 1963 0.875
# 4 1964 3.49 
# 5 1965 3.05 
# 6 1966 4.15 
# 7 1967 4.92 
# 8 1968 11.4 
# 9 1969 9.74 
# 10 1970 8.77 
# # ... with 47 more rows

Visualização de dados

Visualização é uma etapa essencial para entender os dados, o que permite identificar padrões e modelos apropriados. No nosso exemplo, criamos um gráfico de linha para plotar a série do PIB brasileiro usando a função autoplot():

pib_br %>% 
fabletools::autoplot(Growth) +
ggplot2::labs(title = "Crescimento anual do PIB do Brasil", y = "%")

Podemos também plotar os correlogramas ACF e PACF para identificar o processo estocástico da série, obtendo alguns modelos candidatos:

pib_br %>% 
feasts::gg_tsdisplay(Growth, plot_type = "partial")

Especificação do modelo

Existem muitos modelos de séries temporais diferentes que podem ser usados para previsão, e especificar um modelo apropriado para os dados é essencial para produzir previsões.

Os modelos no framework do fable são especificados usando funções com nomenclatura abreviada do nome do modelo (por exemplo, ARIMA(), AR(), VAR(), etc.), cada uma usando uma interface de fórmula (y ~ x). As variáveis de resposta são especificadas à esquerda da fórmula e a estrutura do modelo é escrita à direita.

Por exemplo, um modelo ARIMA(1,0,2) para a taxa de crescimento do PIB pode ser especificado com: ARIMA(Growth ~ pdq(1, 0, 2)).

Neste caso, a variável resposta é Growth e está sendo modelada usando a estrutura de um modelo ARMA(1, 2) especificada na função especial pdq().

Existem diversas funções especiais para definir a estrutura do modelo e em ambos os lados da fórmula pode ser aplicado transformações. Consulte detalhes da documentação do fable.

Estimar o modelo

Identificado um modelo (ou mais) apropriado, podemos em seguida fazer a estimação usando a função model()1.Neste exemplo, estimaremos os seguintes modelos: ARIMA(1,0,2), ARIMA(1,0,0), ARIMA(0,0,2), o algoritmo de seleção automatizada do auto ARIMA criado pelo prof. Rob Hyndman e um passeio aleatório.


fit <- pib_br %>%
fabletools::model(
arima102 = fable::ARIMA(Growth ~ pdq(1, 0, 2)),
arima100 = fable::ARIMA(Growth ~ pdq(1, 0, 0)),
arima002 = fable::ARIMA(Growth ~ pdq(0, 0, 2)),
auto_arima = fable::ARIMA(Growth),
random_walk = fable::RW(Growth)
)

Diagnóstico do modelo

O objeto resultante é uma "tabela de modelo" ou mable, com a saída de cada modelo em cada coluna:


fit

# # A mable: 1 x 5
# arima102 arima100 arima002
# <model> <model> <model>
# 1 <ARIMA(1,0,2) w/ mean> <ARIMA(1,0,0) w/ mean> <ARIMA(0,0,2) w/ mean>
# # ... with 2 more variables: auto_arima <model>, random_walk <model>

Para obter os critérios de informação use a função glance():


fabletools::glance(fit) %>% dplyr::arrange(AICc)

# # A tibble: 5 x 8
# .model sigma2 log_lik AIC AICc BIC ar_roots ma_roots
# <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <list> <list>
# 1 auto_arima 12.8 -150. 307. 307. 313. <cpl [1]> <cpl [1]>
# 2 arima100 12.6 -152. 311. 311. 317. <cpl [1]> <cpl [0]>
# 3 arima102 12.8 -152. 313. 314. 323. <cpl [1]> <cpl [2]>
# 4 arima002 13.1 -153. 314. 315. 322. <cpl [0]> <cpl [2]>
# 5 random_walk 16.0 NA NA NA NA <NULL> <NULL>

Os critérios de informação indicam que, dos modelos estimados, o modelo automatizado ARIMA(1,1,1) apresentou o menor valor de AICc - seguido pelos demais identificados pelos correlogramas ACF e PACF. Com a função gg_tsresiduals() podemos verificar o comportamento dos resíduos deste modelo, indicando que os resíduos se comportam como ruído branco:


fit %>%
dplyr::select(auto_arima) %>%
feasts::gg_tsresiduals()

Um teste de autocorrelação (Ljung Box) retorna um p-valor grande, também indicando que os resíduos são ruído branco:


fabletools::augment(fit) %>%
dplyr::filter(.model == "auto_arima") %>%
fabletools::features(.innov, feasts::ljung_box, lag = 10, dof = 3)

# # A tibble: 1 x 3
# .model lb_stat lb_pvalue
# <chr> <dbl> <dbl>
# 1 auto_arima 8.63 0.281

Também pode ser interessante visualizar o ajuste do modelo. Utilize a função augment() para obter os valores estimados:


fit %>%
fabletools::augment() %>%
dplyr::filter(.model == "auto_arima") %>%
ggplot2::ggplot(ggplot2::aes(x = Year)) +
ggplot2::geom_line(ggplot2::aes(y = Growth, colour = "Observado")) +
ggplot2::geom_line(ggplot2::aes(y = .fitted, colour = "Modelo")) +
ggplot2::scale_colour_manual(
values = c(Observado = "#282f6b", Modelo = "#b22200")
) +
ggplot2::labs(title = "Crescimento anual do PIB do Brasil", colour = NULL)

Previsão

Com o modelo escolhido, previsões podem ser geradas com a função forecast() indicando um horizonte de escolha.


fit %>%
dplyr::select(auto_arima, random_walk) %>%
fabletools::forecast(h = 5) %>%
fabletools::autoplot(pib_br) +
ggplot2::facet_wrap(~.model)

Perceba que os pontos de previsão médios gerados são bastante similares a um processo de passeio aleatório (equivalente a um ARIMA(0,1,0)). O trabalho adicional de especificar termos AR e MA trouxe pouca diferença para os pontos de previsão neste exemplo, apesar de ser perceptível que os intervalos de confiança do modelo auto ARIMA são mais estreitos do que de um passeio aleatório.

Além disso, a previsão fora da amostra gerada ficou bastante aquém dos dados reais para a taxa de crescimento do PIB brasileiro observados no horizonte em questão, configurando apenas um exercício didático.

Saiba mais

Estes são apenas alguns dos recursos e ferramentas disponíveis na família de pacotes do tidyverts. Para uma referência aprofundada, confira o livro Forecasting: Principles and Practice, 3rd Edition, de Hyndman e Athanasopoulos (2021).

Confira outros exercícios aplicados com pacotes do tidyverts:


[1] A função suporta estimação dos modelos com computação paralela usando o pacote future, veja detalhes na documentação e este post para saber mais sobre o tema.

 

Como reverter a primeira diferença de uma série temporal?

By | Data Science

 

Modelagem de séries temporais frequentemente exige a aplicação de transformações nas variáveis, tal como a bem conhecida primeira diferença. Formalmente, podemos descrever essa transformação como:

Ou seja, dado uma série temporal regularmente espaçada, subtraímos do valor em t o valor anterior (t-1), obtendo a série dita "na primeira diferença" ou "nas diferenças".

A mudança de nível da série geralmente contorna diversas características "não desejadas" pelo analista como tendência, sazonalidade, etc., mas dificulta a interpretação dos valores - especialmente quando pretende-se comunicá-los para públicos não técnicos. Para este objetivo é de grande utilidade saber como diferenciar uma série temporal e, sobretudo, também como reverter a transformação.

Para reverter a primeira diferença de uma série utilizamos a soma cumulativa, representada como:

Ou seja, no final teremos um vetor com a soma de cada valor em t com todos os anteriores, valores estes que neste caso podem representar a série na primeira diferença que se deseja reverter.

Exemplo no R

Para exemplificar, vamos aplicar um exercício simples no R com o objetivo de 1) tomar a primeira de uma série e 2) reverter a transformação para obter a série original.

Para esse exemplo você precisará dos seguintes pacotes:


library(magrittr) # CRAN v2.0.1
library(GetBCBData) # CRAN v0.6
library(dplyr) # CRAN v1.0.7
library(timetk) # CRAN v2.6.2
library(tsibble) # CRAN v1.0.1
library(tidyr) # CRAN v1.1.4

Utilizaremos uma série conhecidamente não estacionária: a taxa de câmbio (R$/US$). Primeiro, importamos a série diretamente do banco de dados do Banco Central (SGS/BCB) e tratamos os dados para obter um tibble:


dados <- GetBCBData::gbcbd_get_series(
id = 3696,
first.date = "2000-01-01",
use.memoise = FALSE
) %>%
dplyr::select("date" = "ref.date", "value") %>%
dplyr::as_tibble()
dados

Uma rápida visualização dos dados:


timetk::plot_time_series(
.data = dados,
.date_var = date,
.value = value,
.smooth = FALSE,
.title = "Taxa de câmbio (R$/US$)"
)

Para criar a série na primeira diferença usamos a função tsibble::difference1, com seus argumentos padrão:

dados %<>%
dplyr::mutate(value_diff = tsibble::difference(value))
dados

Quando a diferença é aplicada inevitavelmente perdemos uma observação (a primeira), sendo preenchida no objeto como um NA. Sendo assim, a reversão dessa série para o seu nível (valores originais) requer uma artimanha: devemos substituir o valor NA pelo valor original da série correspondente a esse período dos dados e, então, usar a função cumsum para calcular a soma cumulativa.
dados %<>%
dplyr::mutate(
value_diff = dplyr::if_else(is.na(value_diff), value, value_diff),
value_revert = cumsum(value_diff)
)
dados
Checando se a série revertida corresponde aos valores originais:
all.equal(dados$value_revert, dados$value)
# [1] TRUE

Os cálculos ocorreram conforme o esperado: a série original da taxa de câmbio importada do BCB é igual à série que aplicamos e revertemos a primeira diferença.

Simples, não?

Por fim, vale comparar visualmente o comportamento da série em nível e na diferença (primeiro gráfico):


dados %>%
dplyr::select(
"date",
"Nível" = "value",
"1ª diferença" = "value_diff",
"Série revertida" = "value_revert"
) %>%
tidyr::pivot_longer(
cols = -"date",
names_to = "variable"
) %>%
timetk::plot_time_series(
.date_var = date,
.value = value,
.facet_vars = variable,
.smooth = FALSE,
.title = "Taxa de câmbio (R$/US$)",
.interactive = FALSE
)

 


[1] Para um equivalente em base R veja ?diff.

Gerando previsões desagregadas de séries temporais

By | Data Science

Frequentemente, séries temporais podem ser desagregadas em vários atributos de interesse. Por exemplo1, o número total de veículos produzidos no país ao longo do tempo pode ser dividido em veículos comerciais leves, caminhões, ônibus, etc. Cada uma dessas categorias pode ainda ser desagregada em outros níveis, por regiões, por fabricante, etc., caracterizando o que pode ser chamado de "séries temporais hierárquicas".

Essa riqueza de informação e dados possibilita (e pode ser de interesse) gerar previsões desagregadas das séries, de modo que os pontos de previsão das séries desagregadas possam ser analisados individualmente e que, quando agregados de alguma forma, sejam coerentes com os valores agregados da série.

De modo geral, podemos representar o exemplo como:

Ou seja, para qualquer observação no tempo t, as observações no nível inferior da hierarquia serão agregadas para formar as observações da série acima.

Previsão hierárquica

Existem várias abordagens para gerar previsões hierárquicas com séries temporais, a mais comum e intuitiva é a bottom-up (de baixo para cima), onde primeiro geramos as previsões para cada série no nível inferior e, em seguida, somamos os pontos de previsão para obter previsões da série acima.

Por outro lado, a abordagem top-down (de cima para baixo) envolve primeiro a geração de previsões para a série agregada e, em seguida, desagregá-la na hierarquia utilizando proporções (geralmente baseadas na série histórica).

Ainda existem outras abordagens como a do Minimum Trace e suas variantes, sendo que cada uma tem seus prós e contras. Para se aprofundar no tema veja Forecasting: principles and practice de Hyndman, R.J., & Athanasopoulos, G. (2021), que apresenta detalhes sobre as abordagens disponíveis.

Exemplo: dados de produção da ANFAVEA

Mensalmente a ANFAVEA disponibiliza séries temporais da produção, licenciamento, exportação, etc. de veículos. Os dados são desagregados pelas categorias citadas acima e neste exemplo iremos explorar alguns métodos de gerar modelos de previsão para as séries de produção de veículos.

Primeiro, para reproduzir o código certifique-se de que tenha os seguintes pacotes em sua instalação do R:

# Carregar pacotes/dependências
library(magrittr)
library(GetBCBData)
library(dplyr)
library(tsibble)
library(fabletools)
library(ggplot2)
library(lubridate)
library(fable)

Para importar as séries podemos usar o banco de dados do Banco Central (BCB) através do pacote GetBCBData, bastando apontar os códigos de coleta obtidos diretamente no site da instituição.

raw_anfavea <- GetBCBData::gbcbd_get_series(
id = c(
"Total" = 1373,
"Comerciais leves" = 1374,
"Caminhões" = 1375,
"Ônibus" = 1376
),
first.date = "1993-01-01",
use.memoise = FALSE
)

Além disso fazemos o tratamento necessário para obter um objeto de classe tsibble, pois precisamos desse formato para a próxima etapa de modelagem e previsão:

anfavea <- raw_anfavea %>%
dplyr::select("date" = "ref.date", "variable" = "series.name", "value") %>%
dplyr::mutate(date = tsibble::yearmonth(.data$date)) %>%
tsibble::as_tsibble(key = "variable", index = "date")
anfavea

Uma rápida visualização dos dados sempre é importante:

anfavea %>%
fabletools::autoplot(value) +
ggplot2::facet_wrap(
~variable, 
scales = "free_y"
) +
ggplot2::labs(
title = "Produção de veículos desagregada - Brasil",
x = NULL,
y = "Unidades",
caption = "Dados: ANFAVEA | Elaboração: analisemacro.com.br"
) +
ggplot2::theme(legend.position = "none")

Para gerar previsões hierárquicas usaremos a família de pacotes do tidyverts (você não leu errado), que possibilita um fluxo de trabalho em apenas 5 passos:

1) Transforme os dados para um objeto tsibble contendo as séries desagregadas do nível inferior (em nosso exemplo todas exceto o "Total");

2) Utilize fabletools::aggregate_key() para definir a estrutura de agregação;

3) Especifique um modelo a ser estimado com fabletools::model() para cada série e suas desagregações;

4) Utilize fabletools::reconcile() para indicar uma abordagem (bottom-up, top-down, etc.) para gerar previsões desagregadas coerentes de acordo com o modelo;

5) Por fim, use a função fabletools::forecast() para gerar previsões para toda a estrutura de agregação.

Especificando a hierarquia da série:


prod_veiculos <- anfavea %>%
dplyr::filter(!variable == "Total") %>%
fabletools::aggregate_key(variable, value = sum(value, na.rm = TRUE))
prod_veiculos

Estimar modelos restringindo a amostra e apontar abordagens de previsão hierárquica:


fit <- prod_veiculos %>%
dplyr::filter(lubridate::year(.data$date) <= 2017) %>%
fabletools::model(baseline = fable::ETS(value)) %>%
fabletools::reconcile(
bottom_up = fabletools::bottom_up(baseline),
top_down = fabletools::top_down(baseline),
ols = fabletools::min_trace(baseline, method = "ols")
)
fit

Gerar previsões fora da amostra:


fcst <- fabletools::forecast(fit, h = "2 years")
fcst

Visualizar previsões geradas (produção agregada e desagregada):


fcst %>%
fabletools::autoplot(
dplyr::filter(
prod_veiculos,
dplyr::between(
lubridate::year(.data$date),
2015,
2019
)
),
level = NULL,
size = 1.1
) +
ggplot2::facet_wrap(
~variable,
scales = "free_y"
) +
ggplot2::labs(
title = "Previsão hierárquica da produção de veículos",
y = "Unidades",
x = NULL
)

Calcular métricas de acurácia:


fcst %>%
dplyr::filter(fabletools::is_aggregated(variable)) %>%
fabletools::accuracy(prod_veiculos) %>%
dplyr::arrange(RMSE)

O objeto fcst de classe fable contém as previsões baseline bem como as previsões coerentes de cada abordagem, para as três categorias de veículos desagregadas e para a produção total. Vale enfatizar que as previsões agregadas baseline serão diferentes das previsões agregadas pela abordagem bottom_up, por exemplo.

Dado a escala das séries em diferentes níveis de agregação serem bastante diferentes, devido à agregação, deve-se tomar cuidado ao avaliar as métricas de acurácia. Acima exibimos a acurácia para a série agregada (produção total) conforme cada abordagem/modelo.

Neste caso, agregar as previsões baseline usando a abordagem top-down resultou em previsões mais precisas em comparação com a abordagem bottom-up.

Por fim, vale destacar que o objetivo desse exercício não contempla nenhum aprofundamento ou preocupação com a especificação econométrica, ajuste ou acurácia do(s) modelo(s), focando apenas na implementação prática de previsões hierárquicas.

 

Referências

Hyndman, R.J., & Athanasopoulos, G. (2021) Forecasting: principles and practice, 3rd edition, OTexts: Melbourne, Australia. OTexts.com/fpp3. Accessed on <2021-12-23>.

Athanasopoulos, G., Ahmed, R. A., & Hyndman, R. J. (2009). Hierarchical forecasts for Australian domestic tourism. International Journal of Forecasting, 25, 146–166.

Gross, C. W., & Sohl, J. E. (1990). Disaggregation methods to expedite product line forecasting. Journal of Forecasting, 9, 233–254.

Hyndman, R. J., Ahmed, R. A., Athanasopoulos, G., & Shang, H. L. (2011). Optimal combination forecasts for hierarchical time series. Computational Statistics and Data Analysis, 55(9), 2579–2589.

Panagiotelis, A., Athanasopoulos, G., Gamakumara, P., & Hyndman, R. J. (2021). Forecast reconciliation: A geometric view with new insights on bias correction. International Journal of Forecasting, 37(1), 343–359.

Wickramasuriya, S., Athanasopoulos, G., & Hyndman, R. (2019). Optimal forecast reconciliation for hierarchical and grouped time series through trace minimization. Journal of the American Statistical Association, 114(526), 804–819.


1 Outros exemplos comuns de séries temporais que podem ser desagregadas são: IPCA, IGP-M, Contas Nacionais, produção industrial, consumo de energia elétrica, etc.

 

Receba diretamente em seu e-mail gratuitamente nossas promoções especiais
e conteúdos exclusivos sobre Análise de Dados!

Assinar Gratuitamente