Analise de Séries Temporais

IMPORTAÇÃO DOS DADOS

EXPLICAÇÃO

tictoc::tic()
dados_br <- readr::read_csv2(file='https://dadosabertos.aneel.gov.br/dataset/5e0fafd2-21b9-4d5b-b622-40438d40aba2/resource/b1bd71e7-d0ad-4214-9053-cbd58e9564a7/download/empreendimento-geracao-distribuida.csv',
                             col_types = readr::cols(MdaPotenciaInstaladaKW = readr::col_number()),
                             locale = readr::locale(decimal_mark = ',', grouping_mark = '.'))


dados_br <-  janitor::clean_names(dados_br)
 
tictoc::toc()

## 508.8 sec elapsed

Selecionando as variáveis de interesse referente ao dados de geração distribuida

library(magrittr)
dados_br <- dados_br %>%
   dplyr::filter(sig_tipo_geracao == "UFV") %>% 
  dplyr::select(dth_atualiza_cadastral_empreend,
                dsc_classe_consumo,
                sig_tipo_consumidor,
                sig_uf,
                cod_municipio_ibge,
                mda_potencia_instalada_kw,
                qtd_uc_recebe_credito)

CRIANDO AS SÉRIES TEMPORAIS

EXPLICAÇÃO

Criando as séries mensais da potência acumulada, número de istalações acumulada e unidades que recebem crédito acumulada, com a abrangência nacional.

library(magrittr)
library(dplyr)
library(lubridate)

# Conversão da coluna de datas
dados_br$dth_atualiza_cadastral_empreend <- as.Date(dados_br$dth_atualiza_cadastral_empreend)

# Processar dados
dados_br_mensal <- dados_br %>% 
  dplyr::select(dth_atualiza_cadastral_empreend,
                mda_potencia_instalada_kw,
                qtd_uc_recebe_credito) %>% 
  dplyr::group_by(dth_atualiza_cadastral_empreend) %>% 
  dplyr::summarise(potencia = sum(mda_potencia_instalada_kw),
                   unidades_recebem_credito = sum(qtd_uc_recebe_credito),
                   count = n()) %>%
  dplyr::ungroup() %>% 
  dplyr::mutate(potencia_acumulada = cumsum(potencia),
                unidades_recebem_credito_acu = cumsum(unidades_recebem_credito),
                numero_de_instalacoes_acumuladas = cumsum(count)) %>%
  dplyr::mutate(ref_data = as.Date(paste0(year(dth_atualiza_cadastral_empreend), '-', month(dth_atualiza_cadastral_empreend), '-01'))) %>% 
  dplyr::group_by(ref_data) %>% 
  dplyr::summarise(potencia_acumulada = last(potencia_acumulada),
                   unidades_recebem_credito_acu = last(unidades_recebem_credito_acu),
                   numero_de_instalacoes_acumuladas = last(numero_de_instalacoes_acumuladas)) %>% 
  dplyr::ungroup() %>% 
  dplyr::select(ref_data,
                potencia_acumulada ,
                unidades_recebem_credito_acu,
                numero_de_instalacoes_acumuladas) %>% 
  dplyr::mutate(log_potencia = ifelse(potencia_acumulada > 0, log(potencia_acumulada), NA),
                log_credito = ifelse(unidades_recebem_credito_acu > 0, log(unidades_recebem_credito_acu), NA),
                log_instalacoes = ifelse(numero_de_instalacoes_acumuladas > 0, log(numero_de_instalacoes_acumuladas), NA))
               
# Visualizar os primeiros registros do resultado
head(dados_br_mensal)

## # A tibble: 6 × 7
##   ref_data   potencia_acumulada unidades_recebem_credit…¹ numero_de_instalacoe…²
##   <date>                  <dbl>                     <dbl>                  <int>
## 1 2009-06-01                8.2                         1                      1
## 2 2009-09-01              187.                         24                     22
## 3 2009-10-01              197.                         25                     23
## 4 2010-03-01              227.                         26                     24
## 5 2010-04-01              228.                         27                     25
## 6 2010-08-01              234.                         29                     26
## # ℹ abbreviated names: ¹unidades_recebem_credito_acu,
## #   ²numero_de_instalacoes_acumuladas
## # ℹ 3 more variables: log_potencia <dbl>, log_credito <dbl>,
## #   log_instalacoes <dbl>

GRÁFICOS DAS SÉRIES

EXPLICAÇÃO

Construindo os gráficos da geração distribuida para o Brasil :

library(magrittr)
library(ggplot2)

grafico_potencia_acumulada_br <-dados_br_mensal %>% 
   ggplot() + 
   geom_line(mapping = aes(x = ref_data,
                           y = potencia_acumulada, 
                           group = 1,
                           text = paste0("Data: ",
                            format(ref_data, "%d-%m-%Y"),
                                         "<br> ","Potência Acumulada: ",
                                         scales::number(potencia_acumulada,
                                                        accuracy =1 ,
                                                        decimal.mark = ',',
                                                        big.mark = '.'))),
             color = "#07020D",
             size = 0.5,
             linetype = 1) +
    geom_point(mapping = aes(x = ref_data,
                           y =potencia_acumulada, 
                           group = 1,
                           text = paste0("Data: ",
                            format(ref_data, "%d-%m-%Y"),
                                         "<br> ","Potência Acumulada: ",
                                         scales::number(potencia_acumulada,
                                                        accuracy =1 ,
                                                        decimal.mark = ',',
                                                        big.mark = '.'))),
             color = "#3E4F5D",
             size = 0.1,
             linetype = 1) +
   scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                     decimal.mark = ',',
                                                     big.mark = '.')) +
   labs(title = 'Serie Temporal Da Potência Acumulada - Brasil - Mensal',
        subtitle = 'Número Toral',
        x = 'Data',
        y = 'kWp')

 
plotly::ggplotly(grafico_potencia_acumulada_br,tooltip = "text")

Fazendo o gráfico referente ao Número de Instalações:

library(magrittr)
library(ggplot2)

grafico_numero_istalacoes_acumulada_br<- dados_br_mensal %>% 
   ggplot() + 
   geom_line(mapping = aes(x = ref_data,
                           y = numero_de_instalacoes_acumuladas, 
                           group = 1,
                           text = paste0("Data: ",
                            format(ref_data, "%d-%m-%Y"),
                                         "<br> ","Número de de Intalaçoes Acumulada: ",
                                         scales::number(numero_de_instalacoes_acumuladas,
                                                        accuracy =1 ,
                                                        decimal.mark = ',',
                                                        big.mark = '.'))),
             color = "#07020D",
             size = 0.5,
             linetype = 1) +
    geom_point(mapping = aes(x = ref_data,
                           y =numero_de_instalacoes_acumuladas, 
                           group = 1,
                           text = paste0("Data: ",
                            format(ref_data, "%d-%m-%Y"),
                                         "<br> ","Número de de Instalações Acumulada: ",
                                         scales::number(numero_de_instalacoes_acumuladas,
                                                        accuracy =1 ,
                                                        decimal.mark = ',',
                                                        big.mark = '.'))),
             color = "#3E4F5D",
             size = 0.1,
             linetype = 1) +
   scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                     decimal.mark = ',',
                                                     big.mark = '.')) +
   labs(title = 'Série Do Número de de Instalações Acumulada - Brasil- Mensal',
        subtitle = 'Número Toral',
        x = 'Data',
        y = 'Número Total')

 
plotly::ggplotly(grafico_numero_istalacoes_acumulada_br,tooltip = "text")

Fazendo o gráfico referente ao Número de Unidades que Recebem Crédito Mensal

library(magrittr)
library(ggplot2)

 grafico_unidades_recebem_credito_acumulada_br<- dados_br_mensal %>% 
   ggplot() + 
   geom_line(mapping = aes(x = ref_data,
                           y = unidades_recebem_credito_acu, 
                           group = 1,
                           text = paste0("Data: ",
                            format(ref_data, "%d-%m-%Y"),
                                         "<br> ","Número de Unidades que Recebem Crédito Acumulada: ",
                                         scales::number(unidades_recebem_credito_acu,
                                                        accuracy =1 ,
                                                        decimal.mark = ',',
                                                        big.mark = '.'))),
             color = "#07020D",
             size = 0.5,
             linetype = 1) +
    geom_point(mapping = aes(x = ref_data,
                           y = unidades_recebem_credito_acu, 
                           group = 1,
                           text = paste0("Data: ",
                            format(ref_data, "%d-%m-%Y"),
                                         "<br> ","Número de Unidades que Recebem Crédito Acumulada: ",
                                         scales::number(unidades_recebem_credito_acu,
                                                        accuracy =1 ,
                                                        decimal.mark = ',',
                                                        big.mark = '.'))),
             color = "#3E4F5D",
             size = 0.1,
             linetype = 1) +
   scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                     decimal.mark = ',',
                                                     big.mark = '.')) +
   labs(title = 'Série Da Quantidade de Unidades que Recebem Cédito Acumulada - Brasil- Mensal',
        subtitle = 'Número Toral',
        x = 'Data',
        y = 'Número Total')

 
plotly::ggplotly( grafico_unidades_recebem_credito_acumulada_br,tooltip = "text")

Trasforamando os dados em séries temporais, iremos trabalhar com os dados `a partir de 2015, pois é quando temos novos dados todos os meses.

# selecionando  os dados a partir de 2015
dados_br_mensal <- dados_br_mensal %>% dplyr::slice(51:n())


dados_br_st<-ts(dados_br_mensal, 
            start = c(2015,1), 
            frequency = 12)



dados_br_st<-ts<-xts::as.xts(dados_br_st)

MODELO DE HOLT

EXPLICAÇÃO

Rodando os modelos de Holt(série apenas com tendência) para as princiapais variáveis

library(forecast)

# Para a variável Potência Acumulada :
holt_potencia <- HoltWinters(dados_br_st$potencia_acumulada, gamma = FALSE)
plot(holt_potencia)

# Gráfico das previsões:
fcast_holt_potencia <- forecast(holt_potencia)
autoplot(x =dados_br_st$ref_data,forecast(holt_potencia, h =12))+
  scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                    decimal.mark = ',',
                                                     big.mark = '.')) +
  labs(x = "Ano", y = " Potência Instalada Acumulada(Kw/h)", 
       title = "Previsão para a Potência Instalada Acumulada")

#Calculando a acurácia do modelo:
accuracy(fcast_holt_potencia)

##                    ME     RMSE      MAE       MPE     MAPE        MASE
## Training set 9628.987 28815.93 13835.32 0.8650939 2.604853 0.009878903
##                    ACF1
## Training set 0.06044604

# Para a variável Núemro de Instalações Acumulada: 
holt_instalacoes <- HoltWinters(dados_br_st$numero_de_instalacoes_acumuladas, gamma = FALSE)
plot(holt_instalacoes)

# Gráfico das previsões:
fcast_holt_instalacoes  <- forecast(holt_instalacoes )
autoplot(x =dados_br_st$ref_data,forecast(holt_instalacoes , h =12))+
  scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                    decimal.mark = ',',
                                                     big.mark = '.')) +
  labs(x = "Ano", y = " Número de Instalaçoes Acumulada", 
       title = "Previsão para o Número de Instalaçoes Acumulada")

#Calculando a acurácia do modelo:
accuracy(fcast_holt_instalacoes )

##                   ME     RMSE      MAE       MPE     MAPE       MASE
## Training set 752.539 2957.284 1380.443 0.8923195 1.671878 0.01091528
##                     ACF1
## Training set -0.02723011

# Para a variável Unidades que Recebem Crédito Acumulada:
holt_credito <- HoltWinters(dados_br_st$unidades_recebem_credito_acu , gamma = FALSE)
plot(holt_credito)

# Gráfico das previsões:
fcast_holt_credito <- forecast(holt_credito)
autoplot(x =dados_br_st$ref_data,forecast(holt_credito, h =12))+
  scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                    decimal.mark = ',',
                                                     big.mark = '.')) +
  labs(x = "Ano", y = "Número de Unidades que Recebem Crédito Acumulada", 
       title = "Previsão para o Número Unidades que Recebem Crédito")

#Calculando a acurácia do modelo:
accuracy(fcast_holt_credito)

##                    ME     RMSE      MAE       MPE     MAPE       MASE
## Training set 1150.202 5299.038 2552.741 0.9612625 2.084086 0.01360287
##                     ACF1
## Training set -0.04639439

MODELO DE RNA - FEED-FORWARD

EXPLICAÇÃO

Testando um modelo de previsão RNA para as variáveis:

library(forecast)
#Para a variável Potência Instaltada:

fit_rna_potencia <- nnetar(dados_br_st$potencia_acumulada)
fcast_rna_potencia <- forecast(fit_rna_potencia) 
autoplot(forecast(fit_rna_potencia, h =12)) +
  scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                    decimal.mark = ',',
                                                    big.mark = '.')) +
  labs(x = "Ano", y = "Potência Instalada Acumulada(Kw/h)", 
       title = "Previsão para a Potência Instalada Acumulada")

#Calculando a acurácia do modelo:
accuracy(fcast_rna_potencia)

##                    ME     RMSE      MAE      MPE     MAPE        MASE
## Training set -133.286 18237.94 10784.85 -9.31372 9.777396 0.007700758
##                     ACF1
## Training set -0.02696085

# Para a variável Núemro de Instalações Acumulada: 

fit_rna_instalacoes <- nnetar(dados_br_st$numero_de_instalacoes_acumuladas)
fcast_rna_instalacoes <- forecast(fit_rna_instalacoes)
autoplot(forecast(fit_rna_potencia, h =12)) +
scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                    decimal.mark = ',',
                                                     big.mark = '.')) +
  labs(x = "Ano", y = "Número de Instalaçoes Acumulada", 
       title = "Previsão para o Número de Instalaçoes Acumulada")

#Calculando a acurácia do modelo:
accuracy(fcast_rna_instalacoes)

##                    ME     RMSE      MAE       MPE     MAPE       MASE      ACF1
## Training set 7.643359 2147.658 1285.411 -18.30442 19.01892 0.01016386 0.2204762

# Para a variável Unidades que Recebem Crédito Acumulada:

fit_rna_credito <- nnetar(dados_br_st$unidades_recebem_credito_acu )
fcast_rna_credito <- forecast(fit_rna_credito)
autoplot(forecast(fit_rna_credito, h =12)) +
scale_y_continuous(labels = scales::number_format(accuracy = 1,
                                                    decimal.mark = ',',
                                                     big.mark = '.')) +
  labs(x = "Ano", y = "Número de Unidades que Recebem Crédito Acumulada", 
       title = "Previsão para o Número de Unidades que Recebem Crédito Acumulada")

#Calculando a acurácia do modelo:
accuracy(fcast_rna_credito)

##                      ME     RMSE      MAE       MPE     MAPE       MASE
## Training set -0.5148085 3683.755 2092.705 -22.83939 23.78226 0.01115147
##                    ACF1
## Training set 0.06220853

MODELO DE RNA - LSTM

EXPLICAÇÃO

library(TSLSTMplus)