Research Macro - IPCA
Lucca Pavan
12/04/2021
1 Visualização dos dados
inicio <- Sys.Date()-3660
ipcadata <- BETSget(433, from = paste(inicio), to = paste(Sys.Date()), frequency = 12, data.frame = TRUE)
ipcaxts <- xts(ipcadata$value, order.by = ipcadata$date)
colnames(ipcaxts) <- "ipca"
ipcaxts$ipcatrim <- acum_p(ipcaxts$ipca,3)
ipcaxts$ipcaanual <- acum_p(ipcaxts$ipca,12)
ipcaxts$ipcaanual4 <- lag.xts(ipcaxts$ipcaanual, -4)
ipcaxts$y <- ipcaxts$ipcaanual4 - ipcaxts$ipcaanual
ipcaxts$x1 <- ipcaxts$ipcaanual - lag.xts(ipcaxts$ipcaanual, 1)
ipcaxts$x2 <- lag.xts(ipcaxts$x1,1)
ipcaxts$x3 <- lag.xts(ipcaxts$x2,1)
ipcaxts$x4 <- lag.xts(ipcaxts$x3,1)
ipcaxts <- ipcaxts%>%na.omit()plot(ipcaxts[,1:5], legend.loc = "topright", main = "Medidas de inflação \n e variável de interesse, 'y'.")
2 Modelo 1: Recursivo
O primeiro modelo prevê a inflação (IPCA) um ano à frente usando lags do próprio IPCA (mais especificamente, valores passados da mudança percentual trimestral anualizada do IPCA). A estimação é feita na forma recursiva, incluindo 120 observações.
Seguindo (Stock and Watson 1999), o modelo a ser estimado é:
\[\begin{equation} \pi^{4}_{t+4} - \pi_t = \alpha + \beta(L)(\pi_t - \pi_{t-1}) + \epsilon_t \end{equation}\]
Onde, \(\pi^{4}_{t+4}\) é a inflação anual 4 trimestres à frente, \(\pi_t\) é a taxa de inflação trimestral anualizada. \(\beta(L)\) é o operador de lag polinomial.
modelo1 <- lm(y ~x1+x2+x3+x4, data = ipcaxts)
summary(modelo1)##
## Call:
## lm(formula = y ~ x1 + x2 + x3 + x4, data = ipcaxts)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.83065 -0.65883 0.09187 0.80032 1.99837
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0003107 0.1189908 -0.003 0.99792
## x1 0.9328022 0.3131743 2.979 0.00367 **
## x2 -0.1855759 0.3554217 -0.522 0.60278
## x3 0.2066942 0.3568832 0.579 0.56383
## x4 -0.0010989 0.3198012 -0.003 0.99727
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.193 on 96 degrees of freedom
## Multiple R-squared: 0.1013, Adjusted R-squared: 0.0639
## F-statistic: 2.707 on 4 and 96 DF, p-value: 0.03473newdata <- xts(ipcadata$value, order.by = ipcadata$date)
colnames(newdata) <- "ipca"
newdata$ipcaanual <- acum_p(newdata$ipca,12)
newdata$x1 <- newdata$ipcaanual - lag.xts(newdata$ipcaanual,1)
predata <- tail(newdata$x1, 6)
pred.modelo1 <- predict.lm(lm(y ~x1, data = ipcaxts), newdata = predata, se.fit = TRUE)
pred.modelo1## $fit
## 2020-10-01 2020-11-01 2020-12-01 2021-01-01 2021-02-01 2021-03-01
## 0.67742004 0.33853096 0.17642098 0.03351055 0.54988936 0.78252186
##
## $se.fit
## 2020-10-01 2020-11-01 2020-12-01 2021-01-01 2021-02-01 2021-03-01
## 0.2409210 0.1582991 0.1301597 0.1178888 0.2075756 0.2695820
##
## $df
## [1] 99
##
## $residual.scale
## [1] 1.177844accuracy(modelo1)## ME RMSE MAE MPE MAPE MASE
## Training set -1.68996e-17 1.163138 0.9273248 -485.9158 696.2241 0.9258647coeftest(modelo1)##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.00031072 0.11899078 -0.0026 0.997922
## x1 0.93280216 0.31317434 2.9785 0.003667 **
## x2 -0.18557595 0.35542167 -0.5221 0.602783
## x3 0.20669423 0.35688321 0.5792 0.563834
## x4 -0.00109891 0.31980123 -0.0034 0.997265
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Previsão usando pacote forecast
arima_ipca <- auto.arima(ipcaxts$y, xreg = ipcaxts[,6:9])
summary(arima_ipca)## Series: ipcaxts$y
## Regression with ARIMA(2,0,1) errors
##
## Coefficients:
## ar1 ar2 ma1 x1 x2 x3 x4
## 1.3405 -0.4432 0.2865 -0.6616 -0.5635 -0.2849 -0.0199
## s.e. 0.1569 0.1570 0.1795 0.1064 0.1449 0.1500 0.1127
##
## sigma^2 estimated as 0.1624: log likelihood=-49.5
## AIC=115.01 AICc=116.57 BIC=135.93
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.008619768 0.388794 0.2945157 -154.2828 261.7365 0.6803726
## ACF1
## Training set 0.007585775checkresiduals(arima_ipca)
##
## Ljung-Box test
##
## data: Residuals from Regression with ARIMA(2,0,1) errors
## Q* = 9.66, df = 3, p-value = 0.02169
##
## Model df: 7. Total lags used: 103 Modelo 2: Ingênuo
A previsão da inflação um ano à frente é simplesmente a taxa de crescimento anualizada dos mesmos quatro trimestres no ano anterior
\[\begin{equation} \pi^4_{t+4} = \pi^4_{t} = \Big(\prod_{i=1}^4 (1+\pi_{t-i})\Big)^{1/4}-1 \end{equation}\]
modelo2 <- tibble(Data = time(tail(ipcaxts$ipcaanual))+120,`Previsão Ingênua` = tail(ipcaxts$ipcaanual))
print(modelo2)## # A tibble: 6 x 2
## Data `Previsão Ingênua`[,"ipcaanual"]
## <date> <xts>
## 1 2020-09-29 2.132156
## 2 2020-10-29 2.305451
## 3 2020-11-29 2.438302
## 4 2020-12-30 3.135162
## 5 2021-01-29 3.918206
## 6 2021-03-01 4.311091RMSE_ingenuo <- sum((ipcaxts$ipcaanual4 - ipcaxts$ipcaanual)^2)^(1/2)
print(RMSE_ingenuo)## [1] 12.331164 Modelo 3 : Curva de Phillips
Consiste no Modelo 1 adicionando medidas de atividade econômica e lags de inflação escolhendo o " lag ótimo" por meio do Critério Bayesiano de Informação (BIC). a equação de regressão é:
\[\begin{equation} \pi^{4}_{t+4} - \pi_t = \alpha + \beta(L)(\pi_t - \pi_{t-1}) + \gamma(L) x_t + \epsilon_t \end{equation}\]
Em que , \(x_t\) reprsentam as medidas de atividade econômica. Neste caso, usarei o hiato PIB como variável de atividade econômica. Como proxy para o PIB mensal usarei o indice IBC-Br do Banco Central. Para calcular o hiato usarei o resíduo de um modelo ARIMA da série.
ibc <- BETSget(24364, from = paste(inicio), to = paste(Sys.Date()), frequency = 12, data.frame = TRUE)
ibcxts <- xts(ibc[,2], order.by = ibc[,1], frequency = 12)
colnames(ibcxts) <- "ibc"
autoplot(ibcxts)
ibc_arima <-auto.arima(ibcxts, seasonal = TRUE)
summary(ibc_arima)## Series: ibcxts
## ARIMA(0,1,1)
##
## Coefficients:
## ma1
## 0.3009
## s.e. 0.0862
##
## sigma^2 estimated as 3.248: log likelihood=-236.47
## AIC=476.95 AICc=477.05 BIC=482.49
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0201108 1.786898 1.056648 0.007679281 0.7779474 1.0039
## ACF1
## Training set -0.009230673checkresiduals(ibc_arima)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,1)
## Q* = 6.6791, df = 9, p-value = 0.6705
##
## Model df: 1. Total lags used: 10# hiato do ibc
ibc_hiato <- xts(ibc_arima[["residuals"]], order.by = ibc$date)
autoplot(ibc_hiato)
ibc_hiato <- cbind(ibc_hiato, lag.xts(ibc_hiato, c(1:4)))
colnames(ibc_hiato) <- c("ibc_hiato",paste("ibc_hiato lag", c(1:4)))ipcaxts<- merge.xts(ipcaxts,ibc_hiato)%>%na.omit()
ipca_phillips <- auto.arima(ipcaxts$y, xreg = ipcaxts[,6:ncol(ipcaxts)])
summary(ipca_phillips)## Series: ipcaxts$y
## Regression with ARIMA(2,0,0) errors
##
## Coefficients:
## ar1 ar2 x1 x2 x3 x4 ibc_hiato
## 1.5312 -0.6234 -0.6283 -0.5409 -0.2692 -0.0202 -0.0066
## s.e. 0.0789 0.0831 0.1123 0.1424 0.1477 0.1222 0.0210
## ibc_hiato.lag.1 ibc_hiato.lag.2 ibc_hiato.lag.3 ibc_hiato.lag.4
## -0.0299 -0.0400 -0.0508 -0.0401
## s.e. 0.0321 0.0339 0.0306 0.0186
##
## sigma^2 estimated as 0.1656: log likelihood=-48.27
## AIC=120.55 AICc=124.09 BIC=151.93
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.00548558 0.3841343 0.2957745 -694.5407 788.0532 0.6832806
## ACF1
## Training set 0.07592742coeftest(ipca_phillips)##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 1.5312160 0.0788612 19.4166 < 2.2e-16 ***
## ar2 -0.6233627 0.0830834 -7.5029 6.244e-14 ***
## x1 -0.6283217 0.1123120 -5.5944 2.213e-08 ***
## x2 -0.5408580 0.1424311 -3.7973 0.0001463 ***
## x3 -0.2692060 0.1476612 -1.8231 0.0682832 .
## x4 -0.0202217 0.1222418 -0.1654 0.8686102
## ibc_hiato -0.0065555 0.0210485 -0.3114 0.7554601
## ibc_hiato.lag.1 -0.0298547 0.0320780 -0.9307 0.3520125
## ibc_hiato.lag.2 -0.0399910 0.0339306 -1.1786 0.2385528
## ibc_hiato.lag.3 -0.0508122 0.0306433 -1.6582 0.0972809 .
## ibc_hiato.lag.4 -0.0400579 0.0185726 -2.1568 0.0310190 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ipca_phillips <- Arima(ipcaxts$y, order = c(2,0,0), xreg = ipcaxts[,6:ncol(ipcaxts)])
summary(ipca_phillips)## Series: ipcaxts$y
## Regression with ARIMA(2,0,0) errors
##
## Coefficients:
## ar1 ar2 intercept x1 x2 x3 x4 ibc_hiato
## 1.5304 -0.6222 0.0664 -0.6289 -0.5418 -0.2703 -0.0211 -0.0064
## s.e. 0.0790 0.0834 0.4084 0.1123 0.1424 0.1477 0.1223 0.0210
## ibc_hiato.lag.1 ibc_hiato.lag.2 ibc_hiato.lag.3 ibc_hiato.lag.4
## -0.0297 -0.0398 -0.0507 -0.0400
## s.e. 0.0321 0.0339 0.0306 0.0186
##
## sigma^2 estimated as 0.1674: log likelihood=-48.26
## AIC=122.52 AICc=126.71 BIC=156.52
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -0.0006296728 0.3840865 0.2951292 -657.3437 749.9286 0.68179
## ACF1
## Training set 0.076336595 Modelo 4: Núcleos de inflação
Este modelo usa as diferentes medidas de inflação núcleo como regressores, como indicado em Pasaogullari and Meyer (2010) e fica:
\[\begin{equation} \pi^4_{t+4} = \alpha + \lambda(L)\pi_t^* + \epsilon_t \end{equation}\]
Em que \(\pi^*_t\) é o tipo de medida de inflação núcleo. Foram usados no máximo 4 defasagens e a quantidade ótima de defasagens foi determinada pelo critério Bayesiano de informação.
## Pegar núcleos
codes = c(4466,11426,11427,16121,16122, 27838, 27839)
nucleos = BETSget(codes, from='2006-07-01', data.frame=T)## Warning:
## BETS-package: There is no corresponding entry in the metadata table.
##
## Don't worry, this is not a critical problem. We are working on a solution.
## Warning:
## BETS-package: There is no corresponding entry in the metadata table.
##
## Don't worry, this is not a critical problem. We are working on a solution.data_nucleos = matrix(NA, nrow=nrow(nucleos[[1]]),
ncol=length(codes))
for(i in 1:length(codes)){
data_nucleos[,i] = t(nucleos[[i]]$value)
}
colnames(data_nucleos) = c('ipca_ms', 'ipca_ma', 'ipca_ex0',
'ipca_ex1', 'ipca_dp', 'ipca_ex2',
'ipca_ex3')
library(tidyverse)
library(tstools)
nucleos_vm =
data_nucleos %>%
as_tibble() %>%
mutate(date = nucleos[[1]]$date) %>%
select(date, everything())
nucleos_12m =
data_nucleos %>%
ts(start=c(2006,07), freq=12) %>%
acum_p(12) %>%
as_tibble() %>%
mutate(date = nucleos[[1]]$date) %>%
select(date, everything()) %>%
drop_na()nucleoxts <- xts(nucleos_12m[,-1], order.by = nucleos_12m$date)
nucleodata <- merge.xts(ipcaxts$ipcaanual4, nucleoxts, lag.xts(nucleoxts, c(1:4)))%>%na.omit()
modelo4 <- auto.arima(nucleodata[,1], xreg = nucleodata[,-1])
summary(modelo4)## Series: nucleodata[, 1]
## Regression with ARIMA(2,0,1) errors
##
## Coefficients:
## ar1 ar2 ma1 intercept ipca_ms ipca_ma ipca_ex0 ipca_ex1
## 1.8146 -0.9318 -0.3965 2.7734 0.4888 1.8866 -0.5225 0.0690
## s.e. 0.0523 0.0511 0.1351 0.6778 0.5293 0.5034 0.2279 0.1657
## ipca_dp ipca_ex2 ipca_ex3 ipca_ms.4 ipca_ma.4 ipca_ex0.4 ipca_ex1.4
## -0.8615 -0.2671 -0.2418 1.6726 1.0794 -0.0235 0.1202
## s.e. 0.3982 0.7239 0.6671 0.5176 0.4843 0.2182 0.1736
## ipca_dp.4 ipca_ex2.4 ipca_ex3.4 ipca_ms.1 ipca_ma.1 ipca_ex0.1
## -1.8979 -2.8702 1.3954 -1.0114 1.7974 -0.8660
## s.e. 0.3971 0.7561 0.7199 0.5383 0.4590 0.1996
## ipca_ex1.1 ipca_dp.1 ipca_ex2.1 ipca_ex3.1 ipca_ms.2 ipca_ma.2
## 0.1912 -1.0358 1.1249 0.1366 -0.8131 1.0398
## s.e. 0.1687 0.3777 0.8335 0.8219 0.5248 0.4265
## ipca_ex0.2 ipca_ex1.2 ipca_dp.2 ipca_ex2.2 ipca_ex3.2 ipca_ms.3
## -0.7028 0.5293 -0.8409 -0.0687 1.2498 -1.0729
## s.e. 0.1976 0.1687 0.3620 0.7348 0.7149 0.6062
## ipca_ma.3 ipca_ex0.3 ipca_ex1.3 ipca_dp.3 ipca_ex2.3 ipca_ex3.3
## 1.3612 -0.6225 0.2038 -0.3235 0.6003 0.3207
## s.e. 0.4896 0.1982 0.1669 0.3589 0.7719 0.8080
##
## sigma^2 estimated as 0.09233: log likelihood=-0.69
## AIC=81.38 AICc=136.04 BIC=185.98
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.004205347 0.2380654 0.1926467 -0.3914171 4.184874 0.5668245
## ACF1
## Training set 0.01852905coeftest(modelo4)##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 1.814579 0.052327 34.6779 < 2.2e-16 ***
## ar2 -0.931817 0.051107 -18.2328 < 2.2e-16 ***
## ma1 -0.396492 0.135064 -2.9356 0.0033291 **
## intercept 2.773435 0.677839 4.0916 4.284e-05 ***
## ipca_ms 0.488770 0.529321 0.9234 0.3558040
## ipca_ma 1.886587 0.503388 3.7478 0.0001784 ***
## ipca_ex0 -0.522456 0.227931 -2.2922 0.0218961 *
## ipca_ex1 0.069045 0.165730 0.4166 0.6769631
## ipca_dp -0.861482 0.398248 -2.1632 0.0305275 *
## ipca_ex2 -0.267119 0.723862 -0.3690 0.7121127
## ipca_ex3 -0.241776 0.667076 -0.3624 0.7170224
## ipca_ms.4 1.672621 0.517587 3.2316 0.0012311 **
## ipca_ma.4 1.079431 0.484293 2.2289 0.0258218 *
## ipca_ex0.4 -0.023460 0.218200 -0.1075 0.9143796
## ipca_ex1.4 0.120159 0.173610 0.6921 0.4888626
## ipca_dp.4 -1.897946 0.397055 -4.7801 1.752e-06 ***
## ipca_ex2.4 -2.870217 0.756079 -3.7962 0.0001469 ***
## ipca_ex3.4 1.395404 0.719865 1.9384 0.0525715 .
## ipca_ms.1 -1.011430 0.538325 -1.8788 0.0602654 .
## ipca_ma.1 1.797363 0.458993 3.9159 9.007e-05 ***
## ipca_ex0.1 -0.866002 0.199556 -4.3396 1.427e-05 ***
## ipca_ex1.1 0.191171 0.168710 1.1331 0.2571588
## ipca_dp.1 -1.035797 0.377703 -2.7424 0.0060999 **
## ipca_ex2.1 1.124863 0.833469 1.3496 0.1771391
## ipca_ex3.1 0.136617 0.821851 0.1662 0.8679753
## ipca_ms.2 -0.813101 0.524764 -1.5495 0.1212708
## ipca_ma.2 1.039790 0.426528 2.4378 0.0147769 *
## ipca_ex0.2 -0.702789 0.197572 -3.5571 0.0003749 ***
## ipca_ex1.2 0.529293 0.168680 3.1378 0.0017020 **
## ipca_dp.2 -0.840896 0.362013 -2.3228 0.0201880 *
## ipca_ex2.2 -0.068679 0.734817 -0.0935 0.9255353
## ipca_ex3.2 1.249835 0.714905 1.7483 0.0804203 .
## ipca_ms.3 -1.072932 0.606219 -1.7699 0.0767478 .
## ipca_ma.3 1.361219 0.489630 2.7801 0.0054343 **
## ipca_ex0.3 -0.622508 0.198215 -3.1406 0.0016862 **
## ipca_ex1.3 0.203789 0.166891 1.2211 0.2220496
## ipca_dp.3 -0.323516 0.358900 -0.9014 0.3673705
## ipca_ex2.3 0.600327 0.771910 0.7777 0.4367364
## ipca_ex3.3 0.320670 0.808037 0.3969 0.6914777
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 16 Conclusão
rbind(c(accuracy(modelo1), "-"), accuracy(ipca_phillips), accuracy(modelo4))## ME RMSE MAE
## "-1.68995988133062e-17" "1.1631379659222" "0.927324848801362"
## Training set "-0.000629672829908862" "0.384086485702373" "0.295129232380045"
## Training set "0.00420534726734043" "0.238065351517412" "0.192646718453191"
## MPE MAPE MASE
## "-485.915789873733" "696.224145860261" "0.925864654534697"
## Training set "-657.343736167758" "749.928632470003" "0.681790021070219"
## Training set "-0.39141709284012" "4.18487395180893" "0.56682450688205"
## ACF1
## "-"
## Training set "0.0763365921209692"
## Training set "0.0185290489081704"O modelo de Núcleos de inflação se mostrou com melhor acurácia entre os modelos para o período estudado.
Referências
Pasaogullari, Mehmet, and Brent Meyer. 2010. “Simple Ways to Forecast Inflation: What Works Best?” Economic Commentary, nos. 2010-17 (December). https://www.clevelandfed.org/newsroom-and-events/publications/economic-commentary/economic-commentary-archives/2010-economic-commentaries/ec-201017-simple-ways-to-forecast-inflation-what-works-best.
Stock, James H, and Mark W Watson. 1999. “Forecasting Inflation.” Journal of Monetary Economics 44 (2): 293–335. https://doi.org/10.1016/S0304-3932(99)00027-6.