Trabalho Final - Análise de Séries Temporais
Série IGP-M
Igor Mazzeto Resende Soares
30/06/2022
Índice Geral de Preços - Mercado
Concebido no final dos anos 1940 para ser uma medida abrangente do movimento de preços
Apura informações sobre a variação de preços do dia 21 do mês anterior ao dia 20 do mês de coleta
Divulgado pelo Instituto Brasileiro de Economia da Fundação Getulio Vargas (FGV IBRE)
Composição:
| Índices | Peso |
|---|---|
| Índice de Preços ao Produtor Amplo (IPA) | 60% |
| Índice de Preços ao Consumidor (IPC)1 | 30% |
| Índice Nacional de Custo da Construção (INCC) | 10% |
Série: IGP-M | Publicação: BACEN
IGP_M = rbcb::get_series(c(igpm = 189), last =207)
Mes = month(IGP_M$date)
Ano = year(IGP_M$date)
IGPM = data.frame(Ano,Mes,IGP_M[,2])
p=12 #frequencia
inicial=c(2005,1)
y = IGPM$igpm
y=ts(y,frequency=p,start=inicial)
autoplot(y,
main = "IGP-M: 2005 a 2020",
ylab = "IGP-M (%)",
xlab = "Anos",
) +
labs(caption = "Fonte: FGV IBRE") +
geom_hline(aes(yintercept = mean(y)),
col = "Red",
linetype = "dashed")+
scale_linetype_manual(name = "Legenda", values = c(2, 2),
guide = guide_legend(override.aes = list(color = c("red","blue"))))
- Aparentemente a série oscila em torno da média
- A partir de jun/2020 observa-se tendência de crescimento da série
- 2021 apresenta alta variabilidade na série
Análise exploratória: Box Plots
par(mfrow = c(1, 2))
boxplot(y~Ano,
main="IGP-M: 2005 a 2020 - Anos",
ylab = "IGP-M (%)",
xlab = "Anos")
abline(h=mean(y), col = "Red", lty = 5)
abline(h=median(y), col = "Blue", lty = 3)
legend("topleft", inset=.02, title="Legenda",
c("Média","Mediana"),col = c("Red","Blue"),lty = c(5,3), cex=0.8)
boxplot(y~Mes,
main="IGP-M: 2005 a 2020 - Meses",
ylab = "IGP-M (%)",
xlab = "Meses")
abline(h=mean(y), col = "Red", lty = 5)
abline(h=median(y), col = "Blue", lty = 3)
legend("topleft", inset=.02, title="Legenda",
c("Média","Mediana"),col = c("Red","Blue"),lty = c(5,3), cex=0.8)
- Alta variação em 2020 e 2021 devido à pandemia (crise nas cadeias produtivas)
- Alta variação nos anos de 2008, 2014 e 2017 (crise das subprimes, Gov. Dilma Roussef, Gov. Michel Temer)
- Nota-se outlier para o mês de setembro
Medidas, histograma e correlogramas
# Medidas descritivas
d = as.data.frame(t(t(c(summary(y),Var=var(y),Sd=sd(y),
skewness=skewness(y),kurtose=kurtosis(y)))))
b = data.frame(round(d,2))
colnames(b)[1] = 'IGP-M'
#Histograma
g1 = ggplot(y, aes(y), main = "Histograma")+
geom_histogram(binwidth =2 * IQR(y) / (length(y)^(1/3))) +
labs( x = "series",
y = "Freq.",
title ="Histograma",
subtitle = "IGP-M",
caption = "Fonte: FGV IBRE")
b1 = ggtexttable(b)
# Correlogramas
y1=window(y, end=c(2021,3))
ap1 = autoplot(y1,main="Periodo Amostral")
dif.1 = ndiffs(y1)
acf1 = ggAcf(y1,lag=24, main = "FAC - Per. Amostral")
pacf1 = ggPacf(y1,lag=24, main = "FACP - Per. Amostral")
grid.arrange(b1,g1,acf1,pacf1,ncol = 2, nrow = 2)
- Distribuição com características leptocúrticas e assimetria para direita
- Baixa variância e desvio-padrão
- Evidência de modelo AR(1) pelos correlogramas
Modelo AR(1) - Resíduos
M1="AR(1)"
modelo1 <- Arima(y1,order=c(1,0,0),seasonal=list(order=c(0,0,0),period=p),include.constant=T)
sm1 = summary(modelo1)
cf1 = coeftest(modelo1)
fits1=fitted(modelo1)
residuos=residuals(modelo1)
sresiduos=scale(residuos)
apres1 = autoplot(residuos)
acfres1 = ggAcf(residuos,lag=24,ci=0.99)
pacf1res = ggPacf(residuos,ci=0.99)
qmodelo1 = ggplot(mapping = aes(sample = sresiduos))+
stat_qq_band()+
stat_qq_point(size = 2)+
stat_qq_line()+
xlab("Theoretical quantiles")+
ylab("Sample quantiles")+
ggtitle("Normal Q-Q Plot of residuals")
ggarrange(qmodelo1,acfres1, pacf1res, ncol = 3)
Modelo AR(1) - Ljung Box
FitAR::LBQPlot(residuos)
Modelo AR(1) - Ajuste do modelo - Gráfico
# Previsao Modelo 1
previsao_M1 <- forecast(modelo1, h=p,level=95)
autoplot(y1, main = "Ajuste do modelo",)+autolayer(fits1)+theme(legend.position="bottom")
Modelo AR(1) - Ajuste do modelo - Resultados
Modelo AR(1) - Previsão
# Modelo real
y_real=window(y, start= c(2021,4))
# Periodo Amostral
#forecast::accuracy(previsao_M1$fitted,y1)
# Periodo real
#forecast::accuracy(previsao_M1$mean,y_real)
# Graficos
aj1 = autoplot(y1, main = "Ajuste do modelo",)+autolayer(fits1)+theme(legend.position="bottom")
prev1 = autoplot(previsao_M1)+autolayer(fits1)+autolayer(y) + theme(legend.position="bottom")
prev1_ic = autoplot(previsao_M1$mean,
main=M1,
ylab = "Previsao",
series = "Previsao")+
autolayer(y_real,series="observado")+
autolayer(cbind(previsao_M1$lower,previsao_M1$upper),series=c("IC 95%","IC 95%"),linetype = 'dashed')+
theme(legend.position="bottom")
ggarrange(prev1,prev1_ic,nrow = 2)
Diferenciação - motivação modelo SARIMA
ndif = ndiffs(y1)
dif_y=diff(y1,lag=ndif,dif=1)
auto_dif = autoplot(dif_y,main="Grafico de Linha")
acf_dif = ggAcf(dif_y,lag=24)
pacf_dif = ggPacf(dif_y,lag=24)
ggarrange(auto_dif,ggarrange(acf_dif, pacf_dif,ncol = 2),nrow= 2)
1 diferenciação
Evidência de médias móveis de segunda ordem
Evidência de média móvel sazonal de primeira ordem
Pelos correlogramas podemos aplicar modelo SARIMA(0,1,2)(0,0,1) - truncagem no PACF
Modelo SARIMA(0,1,2)(0,0,1) - Resíduos
Modelo SARIMA(0,1,2)(0,0,1)- Ljung Box
Modelo SARIMA(0,1,2)(0,0,1) - Ajuste do modelo - Gráfico
Modelo SARIMA(0,1,2)(0,0,1) - Ajuste - Resultados
Modelo SARIMA(0,1,2)(0,0,1) - Previsão
Modelo Holt-Winters - Ajuste
Modelo Holt-Winters - Ajuste do modelo - Resultados
Modelo Holt-Winters - Previsão
Resumo dos modelos ajustados
| Estatísticas | AR(1) | SARIMA(0,2,1)(0,01) | Holt - Winters |
|---|---|---|---|
| \(constante\) | 0.590 | ||
| \(\phi_1\) | 0.658 | ||
| \(\theta_1\) | -0.374 | ||
| \(\theta_2\) | -0.312 | ||
| \(\Theta_1\) | -0.185 | ||
| \(\alpha\) | 0.522 | ||
| \(\delta\) | 1e-04 | ||
| \(\gamma\) | 1e-04 | ||
| MSE | 0.354 | 0.347 | 0.365 |
| AIC | 357.70 | 355.33 | 865.68 |
| BIC | 357.83 | 355.54 | 921.32 |
| Resíduos | AR(1) | SARIMA(0,2,1)(0,01) | Holt - Winters |
|---|---|---|---|
| Média | 0.010 | 0.031 | -0.006 |
| Desvio-Padrão | 0.597 | 0.590 | 0.606 |
| Valor-P | 0.017 | 0.060 |
Comparativo - Período Amostral
| Estatísticas | AR(1) | SARIMA(0,2,1)(0,01) | Holt - Winters |
|---|---|---|---|
| ME | 0.0013 | 0.031 | -0.005 |
| RMSE | 0.595 | 0.589 | 0.604 |
| MAE | 0.463 | 0.456 | 0.458 |
| MPE | -Inf | -Inf | -Inf |
| MAPE | -Inf | Inf | Inf |
| ACF1 | -0.004 | 0.006 | 0.168 |
| Theil’s U | 0 | 0 | 0 |
Modelos AR(1) e SARIMA(0,2,1)(0,0,1)[12] com medidas de erro próximas no período amostral.
Modelo Holt-Winters apresenta maiores erros para o período amostral.
Comparativo - Período de validação
| Estatísticas | AR(1) | SARIMA(0,2,1)(0,01) | Holt - Winters |
|---|---|---|---|
| ME | 0.415 | -0.658 | 1.071 |
| RMSE | 1.146 | 1.281 | 1.612 |
| MAE | -218.66 | -735.58 | -1026.25 |
| MPE | 312.98 | 800.69 | 1115.96 |
| MAPE | 0.171 | 0.083 | 0.121 |
| Theil’s U | 0.297 | 1.399 | 1.435 |
Grande parte das medidas de erro para o AR(1) no período de validação são as menores observadas, principalmente o Theil’s U.
Modelo ajustado para Holt-Winter apresentou piores medidas de precisão do modelo
Conclusões
Pelas estatísticas dos erros de previsão, o melhor modelo no período amostral é o SARIMA(0,1,2)(0,0,1)[12].
Pelas estatísticas dos erros de previsão, o melhor modelo no período de validação é o AR(1).
Nota-se que a primeira observação da previsão em todos os modelos fica fora do IC de 95%. Devido a alta variabilidade provocada pela pandemia de COVID-19, é um desafio realizar ajustes e previsões bem comportadas que fiquem dentro dos limites estabelecidos e será uma situação que os estatísticos e analistas terão que lidar nas séries temporais.
Pelo princípio da parcimônia, utilizaria para prever o IGP-M, o modelo AR(1). O modelo SARIMA(0,2,1)(0,0,1)[12] possui bom ajuste e assertividade, mas é arrojado.
O modelo de Holt-Winters possui erros demasiadamente grandes (comparados aos outros modelos apresentados) quando o utilizamos para realizar previsões, dessa forma não se mostra muito efetivo.
Coletado nas capitais: Rio de Janeiro, São Paulo, Belo Horizonte, Salvador, Recife, Porto Alegre e Brasília.↩︎