Dados

A série temporal plotada na abaixo consiste em 57 overshorts diários consecutivos de um tanque subterrâneo de gasolina em um posto de gasolina no Colorado. Se \(y_t\) é a quantidade medida de combustível no tanque no final do dia \(t\) e \(a_t\) é a quantidade medida vendida menos a quantidade entregue durante o curso do dia \(t\), então o overshort no final do dia \(t\) é definido como \(x_t = y_t - y_{t−1} +a_t\).

Devido ao erro na medição da quantidade atual de combustível no tanque, a quantidade vendida e a quantidade entregue ao posto, vemos \(Y_t\), \(a_t\) e \(X_t\) como valores observados de algum conjunto de variáveis aleatórias \(Y_t\), \(A_t\) e \(X_t\) para \(t = 1,\cdots,57\). Na ausência de qualquer erro de medição e qualquer vazamento no tanque, cada \(x_t\) seria zero.

oshorts = read.table(file = "http://leg.ufpr.br/~lucambio/CE017/20212S/oshorts.tsm")
head(oshorts)
##    V1
## 1  78
## 2 -58
## 3  53
## 4 -65
## 5  13
## 6  -6

Estudo descritivo

par(mfrow=c(1,1), mar=c(3,3,1,2)+.5, mgp=c(1.6,.6,0), pch=19)
X = ts(oshorts$V1)
plot(X, ylab="Galões", xlab="Dias", col="violet")
points(X)
abline(h=0, lty=1, lwd=2)
grid()

Modelos

Como sempre devemos ter uma ideia de qual modelo propôr. Vejamos as funções autocorrelação e autocorrelação parcial o que nos informam, assim como outras correlações.

library(astsa)
lag1.plot(X, 4)

acf2(X)

##      [,1]  [,2]  [,3]  [,4]  [,5] [,6]  [,7] [,8]
## ACF  -0.5  0.12 -0.21  0.08  0.02 0.12 -0.22 0.25
## PACF -0.5 -0.18 -0.32 -0.26 -0.15 0.04 -0.19 0.12

Parece plausível um modelo \(AR(1)\).

Estimação

modelo01 = arima(X, order=c(1,0,0))
modelo01
## 
## Call:
## arima(x = X, order = c(1, 0, 0))
## 
## Coefficients:
##           ar1  intercept
##       -0.5129    -4.4509
## s.e.   0.1141     4.4196
## 
## sigma^2 estimated as 2518:  log likelihood = -304.22,  aic = 614.44
modelo02 = arima(X, order=c(2,0,0))
modelo02
## 
## Call:
## arima(x = X, order = c(2, 0, 0))
## 
## Coefficients:
##           ar1      ar2  intercept
##       -0.6026  -0.1779    -4.6034
## s.e.   0.1308   0.1322     3.7097
## 
## sigma^2 estimated as 2438:  log likelihood = -303.33,  aic = 614.66

PREVISÃO

Aceitamos como modelo adequado \(AR(1)\) (menor AIC). Podemos então proceder a predizer overshorts futuros.

# horizonte de provisão
m = 3
preditos = predict(modelo01, n.ahead = m)
preditos
## $pred
## Time Series:
## Start = 58 
## End = 60 
## Frequency = 1 
## [1]  1.986064 -7.752596 -2.757274
## 
## $se
## Time Series:
## Start = 58 
## End = 60 
## Frequency = 1 
## [1] 50.17688 56.39276 57.91744
par(mfrow = c(1,1), mar=c(4,3,1,1), mgp=c(1.6,.6,0), pch=19)
ts.plot(X, preditos$pred, col=1:2, lwd=2, ylab="Overshorts", xlab="Dias")
U = preditos$pred+preditos$se; U; L = preditos$pred-preditos$se; L
## Time Series:
## Start = 58 
## End = 60 
## Frequency = 1 
## [1] 52.16295 48.64016 55.16017
## Time Series:
## Start = 58 
## End = 60 
## Frequency = 1 
## [1] -48.19082 -64.14536 -60.67472
xx = c(time(U), rev(time(U))); yy = c(L, rev(U))
polygon(xx, yy, border = 8, col = gray(.6, alpha = .2))
lines(preditos$pred, type="p", col=2)
grid()