Métodos de Previsão – MPAgro/FGV
Código utilizado em aula – Modelo ARIMA
Prof. Dr Luciano Rodrigues
07/07/2018
1 – Carregando a base de dados
1.1 – Carregar arquivo .csv
Vamos iniciar com comando setwd para indicar a pasta de trabalho com o arquivo que será utilizado. Posteriormente, utilizaremos o comando read.csv2 para carregar o arquivo arroz.csv. O arquivo .csv possui uma série de preço anual de arroz, em R$/tonelada.
Os dados do arquivo arroz.csv serão armazenados em um data.frame no objeto base. No processo de carga, especificar corretamente os parâmetros header, sep e dec do comando read.csv2
setwd ("c:\\metodos_previsao")
base = read.csv2("arroz.csv", header = TRUE, sep = ";", dec=",")
1.2 – Converter os preços em uma série de tempo
Vamos utilizar o comando ts para definir os preços como uma série temporal. Essa série será gravada no objeto com o nome parroz.
parroz = ts(base$arroz, start=c(1948), end=c(1992), frequency=1)Vamos construir um gráfico com a série para verificar se os dados foram carregados corretamente.
plot(parroz, main="Preço anual do arroz", ylab = "R$/t", xlab = "período")2 – Processo integrado auto-regressivo de médias móveis (ARIMA)
Vamos adotar a abordagem proposta por Box & Jenkins, realizando os seguintes passos:
- Identificação
- Estimação
- Verificação
- Previsão
Antes de iniciar a identificação do modelo, vamos avaliar a necessidade de realizar alguma transformação na série.
plot(diff(parroz))A partir do gráfico da variação da série, é possível verificar que a variância aumenta com o passar do tempo. Para tentar estabilizar a variância da série, vamos realizra a transformação logarítmica. A nova série será armazenada no objeto lparroz (logarítmo natural do preço do arroz).
lparroz = log(parroz)
plot(lparroz)
plot(diff(lparroz))2.1 – Identificação do modelo
Vamos utilizar a função acf2 do pacote astsa para construir as funções de autocorrelação (FAC) e de autocorrelação parcial (FACP) da série lparroz.
require(astsa)
acf2(lparroz, max.lag=15)O gráfico da série e o comportamento da FAC sugerem que a série não é estacionária.
Há, basicamente, quatro maneiras de observar se série em estudo é ou não estacionária: 1. Análise gráfica 2. Comparar a média e a variância em diferentes períodos de tempo 3. Observar a FAC 4. Realizar testes de Raiz Unitária
Nessa primeira aula iremos utilizar a análise visual gráfica e a avaliação da FAC. Na próxima aula vamos utilizar testes de raiz unitária para avaliar a ordem de integração da série.
Para trabalhar com uma série estacionária, vamos aplicar a primeira diferença nos dados do objeto lparroz e construir a FAC e a FACP.
acf2( diff(lparroz) )A análise da FAC e da FACP sugerem os seguintes modelos:
- ARIMA (0,1,1)
- ARIMA (1,1,0)
- ARIMA (1,1,1)
Como os modelos propostos por Box & jenkins são da década de 1970, o esforço computacional para estimar o sistema era muito grande, portanto, a etapa de identificação era fundamental para se ter um modelo adequado. Atualmente, graças aos avanços computacionais, observar a FAC e a FACP é útil para se ter uma ideia inicial dos modelos que serão avaliados. Como veremos mais adiante, vamos selecionar o modelo que minimize os critérios de informação.
2.2 – Estimação dos modelos identificados
Para estimar os modelos identificados, vamos utilizar a função Arima do pacote forecast. Vamos estimar os três modelos com e sem constante (discutiremos esse aspecto nas próximas aulas).
Com relação ao método de estimação dos parâmetros, usaremos o procedimento default da função do R, que adota o método de Máxima Verossimilhança (Maximum Likelihood-ML).
require(forecast)
arima011 = Arima(lparroz, order=c(0,1,1))
arima011_cte = Arima(lparroz, order=c(0,1,1), include.constant = TRUE)
arima110 = Arima(lparroz, order=c(1,1,0))
arima110_cte = Arima(lparroz, order=c(1,1,0), include.constant = TRUE)
arima111 = Arima(lparroz, order=c(1,1,1))
arima111_cte = Arima(lparroz, order=c(1,1,1), include.constant = TRUE)
Os resultados dos modelos estimados foram armazenados nos objetivos arima011,arima011_cte, etc.
Para selecionar o melhor modelo, vamos identificar aquele que minimiza o Critério de Informação de Akaike (AIC) e o Critério Bayesiano ou de Schwarz (BIC).
arima011$aic
arima011_cte$aic
arima110$aic
arima110_cte$aic
arima111$aic
arima111_cte$aic arima011$bic
arima011_cte$bic
arima110$bic
arima110_cte$bic
arima111$bic
arima111_cte$bicAmbos os critérios indicaram o modelo ARIMA(0,1,1) com constante
2.3 – Verificação ou diagnóstico do modelo
Nesta fase as seguintes características serão analisadas:
Significância dos parâmetros estimados
Análise dos resíduos, especialmente avaliar a ausência de autocorrelação linear
2.3.1 – Significância dos parâmetros estimados
O modelo não deve conter parâmetros em excesso (critério da parcimônia). A existência de parâmetros redundantes será avaliada a partir do teste t de Student, cuja estatística é dada pela razão entre o coeficiente estimado e o seu erro-padrão.
Se o valor o coeficiente estimado for pequeno em relção a seu erro-padrão, temos uma indicação de não significância estatístia e, portanto, é provável que exista superespecificação do modelo.
Regra de decisão a ser adotada: se \(|t| > 2\) então \(H_0: \phi_i=0,\ i=1,..,p\) ou \(H_0:\theta_j,\ j=1,..,q\) é rejeitada (ver detalhes na apresentação utilizada em aula).
arima011_cte
-0.4948/0.13942.3.2 – Análise dos resíduos
Os resíduos do modelo estimado são estimativas do ruído branco. Dessa forma, os mesmos devem se comportar aproximadamente como uma ruído branco se o modelo estiver adequadamente espeficiados. Em particular, seus coeficientes de autocorrelação devem ser estatisticamente iguais a zero.
Para verificar se essa condição é atendida, vamos realizar testes individuais e testes conjuntos para os coeficientes de autocorrelação (verificar detalhes dos testes na apresentação utilizada em aula)
acf2(arima011_cte$residuals)
Box.test(arima011_cte$residuals, type=c("Ljung-Box"), lag=2)
Box.test(arima011_cte$residuals, type=c("Ljung-Box"), lag=4)
Box.test(arima011_cte$residuals, type=c("Ljung-Box"), lag=6)
Todos os testes realizados indicam a ausência de autocorrelação dos resíduos.
De maneira mais direta, o diagnóstico dos resíduos também pode ser realizado utilizando a função checkresiduals do pacote forecast.
checkresiduals(arima011_cte)2.4 – Previsão
A última etapa da metodologia de Box & Jenkins consiste em, a partir do modelo estimado, realizar previsões (verificar detalhes no PPT utilizado em aula).
Para realizar as previsões, iremos uitlizar a função forecast do pacote forecast. Como exercício, vamos realizar previsão para três anos (parâmetro h=3), com intervalo de confiança de 95%
previsao = forecast(arima011_cte, h=3, level=0.95)
plot(previsao)Como trabalhamos com a série em ln, precisamos calcular o antilogaritmo dos valores previstos. Uma das formas de executar essa transformação é dada por:
previsao$mean = exp(previsao$mean)
previsao$x = exp(previsao$x)
previsao$upper = exp(previsao$upper)
previsao$lower = exp(previsao$lower)
previsao$fitted = exp(previsao$fitted)
plot(previsao)