import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.holtwinters import ExponentialSmoothing
from prophet import Prophet
import pmdarima as pm
from sklearn.metrics import mean_absolute_error
from statsmodels.tools.eval_measures import rmse
from statsmodels.stats.diagnostic import acorr_ljungbox
# 1. Importar série do BCB (exemplo: índice volume vendas varejo restrito - código 23427)
# serie = sgs.get({'varejo': 23427}, start='2000-01-01')
# serie.index = pd.to_datetime(serie.index)
# serie = serie.asfreq('MS') # frequência mensal
# serie = serie.fillna(method='ffill')
# Criação da série de datas a partir do número de dias desde 1970-01-01
date_values = [10957, 10988, 11017, 11048, 11078, 11109, 11139, 11170, 11201, 11231, 11262, 11292, 11323, 11354, 11382,
11413, 11443, 11474, 11504, 11535, 11566, 11596, 11627, 11657, 11688, 11719, 11747, 11778, 11808, 11839,
11869, 11900, 11931, 11961, 11992, 12022, 12053, 12084, 12112, 12143, 12173, 12204, 12234, 12265, 12296,
12326, 12357, 12387, 12418, 12449, 12478, 12509, 12539, 12570, 12600, 12631, 12662, 12692, 12723, 12753,
12784, 12815, 12843, 12874, 12904, 12935, 12965, 12996, 13027, 13057, 13088, 13118, 13149, 13180, 13208,
13239, 13269, 13300, 13330, 13361, 13392, 13422, 13453, 13483, 13514, 13545, 13573, 13604, 13634, 13665,
13695, 13726, 13757, 13787, 13818, 13848, 13879, 13910, 13939, 13970, 14000, 14031, 14061, 14092, 14123,
14153, 14184, 14214, 14245, 14276, 14304, 14335, 14365, 14396, 14426, 14457, 14488, 14518, 14549, 14579,
14610, 14641, 14669, 14700, 14730, 14761, 14791, 14822, 14853, 14883, 14914, 14944, 14975, 15006, 15034,
15065, 15095, 15126, 15156, 15187, 15218, 15248, 15279, 15309, 15340, 15371, 15400, 15431, 15461, 15492,
15522, 15553, 15584, 15614, 15645, 15675, 15706, 15737, 15765, 15796, 15826, 15857, 15887, 15918, 15949,
15979, 16010, 16040, 16071, 16102, 16130, 16161, 16191, 16222, 16252, 16283, 16314, 16344, 16375, 16405,
16436, 16467, 16495, 16526, 16556, 16587, 16617, 16648, 16679, 16709, 16740, 16770, 16801, 16832, 16861,
16892, 16922, 16953, 16983, 17014, 17045, 17075, 17106, 17136, 17167, 17198, 17226, 17257, 17287, 17318,
17348, 17379, 17410, 17440, 17471, 17501, 17532, 17563, 17591, 17622, 17652, 17683, 17713, 17744, 17775,
17805, 17836, 17866, 17897, 17928, 17956, 17987, 18017, 18048, 18078, 18109, 18140, 18170, 18201, 18231,
18262]
dates = pd.to_datetime(date_values, origin='1970-01-01', unit='D')
# Valores associados
values = [35.2, 35.6, 39.2, 40.5, 41.6, 40.4, 40.8, 38.7, 37.3, 37.6, 35.6, 47.4, 34.2, 32.2, 38, 37.5, 38.8, 35, 38.4,
40.2, 38.2, 39.3, 36, 46.3, 36.4, 34, 39, 37.8, 38.9, 35.3, 37.2, 38.1, 35.6, 38.3, 35.6, 45.7, 32.2, 31.6, 35.2,
36.8, 37.5, 34.8, 38.4, 38.1, 37, 39, 37.3, 49, 35.8, 35.3, 40.2, 41.3, 43.9, 41.6, 45.9, 42.3, 42.2, 44, 41.3,
56.9, 38.5, 38.5, 45.3, 43.6, 46.2, 44.3, 47.5, 46.5, 46.4, 46, 44.1, 61, 42, 40.2, 44.6, 44.5, 47.8, 45.3, 46.5,
48.5, 47.7, 50.2, 49.3, 64.5, 47, 46.8, 51, 50.5, 55, 51.3, 52.8, 55.3, 54.8, 55.6, 55.3, 72.2, 54.5, 52.1, 56.2,
57.2, 60.8, 56.1, 61.8, 61.6, 59.8, 63.3, 57.7, 77.4, 61.4, 51.9, 57.3, 57.9, 61.9, 57.3, 61.1, 61.1, 60.6, 65.6,
63.5, 83.1, 64.1, 60.2, 67.8, 67.1, 72.8, 68.5, 71.1, 69.3, 69.9, 71.1, 67.9, 92.7, 67.5, 64.8, 69.1, 69.4, 79.6,
70.2, 73.8, 72.5, 71.3, 75.6, 74.7, 100.8, 79.5, 75.7, 82.4, 78, 84.8, 83.2, 84.8, 88.5, 86.3, 91.7, 92.8, 111.4,
92.8, 83.7, 92.5, 88.3, 93.9, 88.8, 96, 95.9, 93.2, 98.3, 100.5, 128.8, 97.2, 90.2, 94.3, 94.4, 101.1, 92, 96.4,
98.2, 97.6, 105.8, 103.1, 129.6, 99.6, 87.8, 97, 94.8, 98.6, 93.4, 98.4, 96.4, 92.5, 100.6, 97.2, 124.5, 91.5,
85.1, 91.6, 88.5, 92.2, 87.4, 90.5, 88.1, 85.2, 89.4, 93.4, 116.9, 90.8, 84, 89.7, 86.3, 90, 87.3, 90.8, 93.5,
93.7, 91.4, 93.5, 114.1, 87.8, 81.1, 94.5, 83.2, 89.9, 88.8, 89.3, 93.7, 93.5, 96.3, 101.3, 118.3, 93.8, 85.2,
90, 86.6, 90, 85.2, 91, 94.4, 93.4, 95.9, 102.6, 116, 94.9]
# Criação do DataFrame
varejoms = pd.DataFrame({'date': dates, 'value': values})
# 2. Separar treino e teste
h = 48
train = varejoms.iloc[:-h]
test = varejoms.iloc[-h:]
# 3. ARIMA
model_arima = pm.auto_arima(train['value'],seasonal=True, m=12,stepwise=True,suppress_warnings=True)
fc_arima = model_arima.predict(n_periods=h)
# 4. ETS
model_ets = ExponentialSmoothing(train['value'], seasonal='add', seasonal_periods=12).fit()
fc_ets = model_ets.forecast(h)
# 5. Combinação simples
fc_comb = (pd.Series(fc_arima, index=test.index) +
pd.Series(fc_ets.values, index=test.index)) / 2 #+
# pd.Series(fc_prophet.values, index=test.index)
# 1. Tamanho do forecast
h = 48
futuro = 12
total_forecast = h + futuro
#
# 2. Forecasts para 60 meses
fc_arima = model_arima.predict(n_periods=total_forecast)
fc_ets = model_ets.forecast(total_forecast)
# 3. Criar datas futuras a partir do final do treino
last_date = train['date'].iloc[-1]
future_dates = pd.date_range(start=last_date + pd.offsets.MonthBegin(1), periods=total_forecast, freq='MS')
# 4. Criar DataFrames com essas datas
fc_arima_series = pd.Series(fc_arima, index=future_dates)
fc_ets_series = pd.Series(fc_ets, index=future_dates)
print(fc_arima[:10]) # Verifique as primeiras 10 previsões
## 193 83.413607
## 194 91.167483
## 195 89.302279
## 196 94.181451
## 197 88.139715
## 198 93.038841
## 199 92.357442
## 200 89.499348
## 201 97.129095
## 202 94.631832
## dtype: float64
print(fc_ets[:10]) # Verifique as primeiras 10 previsões
## 193 83.517814
## 194 92.448886
## 195 90.672297
## 196 95.934751
## 197 90.537434
## 198 95.690674
## 199 95.339464
## 200 93.129704
## 201 100.653711
## 202 97.550534
## dtype: float64
last_date = train['date'].iloc[-1]
print(f"Última data de treino: {last_date}")
## Última data de treino: 2016-01-01 00:00:00
# Gerar datas futuras corretamente
future_dates = pd.date_range(start=last_date + pd.offsets.MonthBegin(1), periods=total_forecast, freq='MS')
print(future_dates) # Verifique as datas geradas
## DatetimeIndex(['2016-02-01', '2016-03-01', '2016-04-01', '2016-05-01',
## '2016-06-01', '2016-07-01', '2016-08-01', '2016-09-01',
## '2016-10-01', '2016-11-01', '2016-12-01', '2017-01-01',
## '2017-02-01', '2017-03-01', '2017-04-01', '2017-05-01',
## '2017-06-01', '2017-07-01', '2017-08-01', '2017-09-01',
## '2017-10-01', '2017-11-01', '2017-12-01', '2018-01-01',
## '2018-02-01', '2018-03-01', '2018-04-01', '2018-05-01',
## '2018-06-01', '2018-07-01', '2018-08-01', '2018-09-01',
## '2018-10-01', '2018-11-01', '2018-12-01', '2019-01-01',
## '2019-02-01', '2019-03-01', '2019-04-01', '2019-05-01',
## '2019-06-01', '2019-07-01', '2019-08-01', '2019-09-01',
## '2019-10-01', '2019-11-01', '2019-12-01', '2020-01-01',
## '2020-02-01', '2020-03-01', '2020-04-01', '2020-05-01',
## '2020-06-01', '2020-07-01', '2020-08-01', '2020-09-01',
## '2020-10-01', '2020-11-01', '2020-12-01', '2021-01-01'],
## dtype='datetime64[ns]', freq='MS')
# Verifique o comprimento de cada uma das previsões
print(len(fc_arima), len(fc_ets))
## 60 60
# Verifique se as previsões e as datas têm o mesmo comprimento
assert len(fc_arima) == len(future_dates), "O comprimento das previsões e das datas futuras deve ser o mesmo!"
assert len(fc_ets) == len(future_dates), "O comprimento das previsões e das datas futuras deve ser o mesmo!"
# Criar as séries de previsões associando corretamente ao índice
fc_arima_series = pd.Series(fc_arima.values, index=future_dates)
fc_ets_series = pd.Series(fc_ets.values, index=future_dates)
# Verifique as primeiras previsões
print(fc_arima_series.head())
## 2016-02-01 83.413607
## 2016-03-01 91.167483
## 2016-04-01 89.302279
## 2016-05-01 94.181451
## 2016-06-01 88.139715
## Freq: MS, dtype: float64
print(fc_ets_series.head())
## 2016-02-01 83.517814
## 2016-03-01 92.448886
## 2016-04-01 90.672297
## 2016-05-01 95.934751
## 2016-06-01 90.537434
## Freq: MS, dtype: float64
# Verifique se há NaN nas previsões associadas
print("ARIMA NaN count:", fc_arima_series.isna().sum())
## ARIMA NaN count: 0
print("ETS NaN count:", fc_ets_series.isna().sum())
## ETS NaN count: 0
# Combinação simples das previsões
fc_comb = (fc_arima_series + fc_ets_series) / 2
# Verifique as primeiras previsões combinadas
print(fc_comb.head())
## 2016-02-01 83.465711
## 2016-03-01 91.808184
## 2016-04-01 89.987288
## 2016-05-01 95.058101
## 2016-06-01 89.338574
## Freq: MS, dtype: float64
# Separar as previsões para teste e futuro
fc_comb_test = fc_comb[:h] # Para comparar com o teste real
fc_comb_future = fc_comb[h:] # Para a previsão futura
# Plotar tudo junto
plt.figure(figsize=(12, 6))
# Plotando a série de treino e teste
plt.plot(train['date'], train['value'], label='Treino', color='black', linewidth=2) # cor azul para o treino
plt.plot(test['date'], test['value'], label='Teste real', color='black', linewidth=2) # cor preta para o teste real
# Plotando as previsões com marcadores
plt.plot(fc_arima_series.index, fc_arima_series.values, label='Previsão ARIMA', linestyle='--', marker='o', markersize=5, color='green', linewidth=2)
plt.plot(fc_ets_series.index, fc_ets_series.values, label='Previsão ETS', linestyle='--', marker='x', markersize=5, color='red', linewidth=2)
plt.plot(fc_comb.index, fc_comb.values, label='Forecast combinado', linestyle='--', marker='^', markersize=5, color='blue', linewidth=2)
# Linha vertical indicando o fim do período de teste
plt.axvline(test['date'].iloc[-1], color='grey', linestyle=':', label='Fim do teste')
# Linha vertical indicando o fim do período de treino
plt.axvline(train['date'].iloc[-1], color='grey', linestyle=':', label='Fim do treino')
# Ajustando o título e a legenda
plt.title("Forecast ARIMA, ETS e Combinado", fontsize=16, fontweight='bold')
plt.legend(fontsize=12)
# Ajuste do layout para não cortar os elementos
plt.tight_layout()
# Exibição do gráfico
plt.show()
