Item a)

Considere os dados até fevereiro de 1976. Quais são as previsões para abril e maio de 1976 utilizando médias móveis de 3 meses, 5 meses e 9 meses?

# Implementando o modelo de médias móveis simples

MMS=function(dados,t,r){
  
  # dados:  dados da série temporal 
  # t:      tempo de referência
  # r:      largura do intervalo
  
  mms=sum(dados[I(t-r+1):t])/r
  
  # Como t=14, estão sendo considerados apenas os dados
  # até fevereiro de 1976.
  
}

# Obtendo as previsões
mms_3meses=MMS(vendas,14,3)
mms_5meses=MMS(vendas,14,5)
mms_9meses=MMS(vendas,14,9)

Item b)

E utilizando suavização exponencial simples com α = 0,1; 0,3; 0,7; e 0,9?

# Implementando o modelo SES

SES=function(dados,t,alfa){
  
  k=0
  soma=0
  
  while(k<=I(t-1)){
    # Expressão geral do modelo
    soma=soma+dados[I(t-k)]*((1-alfa)^k)
    k=k+1
  }
  
  # Finalizando o cálculo pós iterações
  ses=alfa*soma+((1-alfa)^t)*dados[1]  # Obs: z0 chapéu = z1
  
}

# Computando as previsões
ses_01=SES(vendas,14,0.1)
ses_03=SES(vendas,14,0.3)
ses_07=SES(vendas,14,0.7)
ses_09=SES(vendas,14,0.9)

Item c)

Dentre todas essas possibilidades, qual é o método que oferece melhores previsões? (compare com os valores verdadeiros somando as diferenças quadráticas)

# Computando as diferenças quadráticas para MMS
SQ_MMS_3=sum((vendas[15:16]-mms_3meses)^2)
SQ_MMS_5=sum((vendas[15:16]-mms_5meses)^2)
SQ_MMS_9=sum((vendas[15:16]-mms_9meses)^2)

# Computando as diferenças quadráticas para SES
SQ_SES_01=sum((vendas[15:16]-ses_01)^2)
SQ_SES_03=sum((vendas[15:16]-ses_03)^2)
SQ_SES_07=sum((vendas[15:16]-ses_07)^2)
SQ_SES_09=sum((vendas[15:16]-ses_09)^2)

Item d)

Utilizando o R, determine o estimador para α e as previsões correspondentes.

library(forecast)

# Ajustando considerando apenas até fev/76
SES_final=ses(vendas[1:14])

summary(SES_final)

Forecast method: Simple exponential smoothing

Model Information:
Simple exponential smoothing 

Call:
 ses(y = vendas[1:14]) 

  Smoothing parameters:
    alpha = 1e-04 

  Initial states:
    l = 56.7148 

  sigma:  33.2016

     AIC     AICc      BIC 
138.8615 141.2615 140.7786 

Error measures:
                     ME     RMSE      MAE       MPE     MAPE      MASE
Training set 0.00250303 30.73874 28.85827 -52.48745 85.15657 0.9018209
                  ACF1
Training set 0.1235235

Forecasts:
   Point Forecast   Lo 80    Hi 80     Lo 95    Hi 95
15       56.71481 14.1652 99.26441 -8.359196 121.7888
16       56.71481 14.1652 99.26441 -8.359196 121.7888
17       56.71481 14.1652 99.26441 -8.359197 121.7888
18       56.71481 14.1652 99.26441 -8.359197 121.7888
19       56.71481 14.1652 99.26441 -8.359197 121.7888
20       56.71481 14.1652 99.26441 -8.359198 121.7888
21       56.71481 14.1652 99.26441 -8.359198 121.7888
22       56.71481 14.1652 99.26441 -8.359198 121.7888
23       56.71481 14.1652 99.26441 -8.359198 121.7888
24       56.71481 14.1652 99.26441 -8.359199 121.7888

Item a)

Faça um boxplot do número de vendas categorizando pelo mês. Observando o gráfico, há indícios de sazonalidade? Justifique.

Observando o gráfico, podemos perceber que o valor mediano em vendas varia ao longo do ano, apresentando altas vendas em meses de novembro a janeiro e também em abril e maio, por exemplo. Também é possível observar que os meses de junho e julho possuem as vendas medianas mais baixas, além de amplitudes reduzidas. Assim sendo, há indícios da existência de sazonalidade.

Item b)

Aplique o método de Holt ou de Holt-Winters (dependendo da resposta do item anterior), estimando os parâmetros no R, considerando os dados até dezembro de 1999.

O método adequado é o de Holt-Winters, visto que há presença de sazonalidade.

# Recuperando apenas os valores até dez/99
bebida_vendas_=bebida_vendas[1:179]

# Definindo a série temporal
bebida_vendas_ts=ts(bebida_vendas_,start=c(1985,2),frequency=12)

# Ajustando o modelo pelo método de Holt-Winters
SEHW=HoltWinters(bebida_vendas_ts,seasonal = "additive")

# Obtendo as estimativas para os parâmetros
SEHW

Holt-Winters exponential smoothing with trend and additive seasonal component.

Call:
HoltWinters(x = bebida_vendas_ts, seasonal = "additive")

Smoothing parameters:
 alpha: 0.4451027
 beta : 0.01339081
 gamma: 0.3826782

Coefficients:
           [,1]
a   123.5803750
b     0.3449261
s1   17.0825573
s2   -2.8511653
s3  -12.8740005
s4    1.2150383
s5   -0.5002987
s6  -23.6073508
s7  -26.9042637
s8   -9.0058172
s9   -3.4752122
s10   2.8330457
s11  12.2183588
s12  15.0566419

Conforme visualizamos na saída do modelo, as estimativas são: \[ \hat{\alpha} = 0.4451027 \\ \hat{\beta} = 0.01339081 \\ \hat{\gamma} = 0.3826782 \\ \]

Item c)

Obtenha as previsões de janeiro até agosto de 2000. Compare com os valores verdadeiros. Você considera que as previsões são adequadas? Justifique.

# Recuperando os dados verdadeiros do ano 2000, de janeiro a agosto
bebida_2000=ts(bebida[180:187,2])

# Recuperando as previsões de janeiro a agosto, do ano 2000
bebida_2000_prev=ts(forecast(SEHW)$mean[1:8])

Vejamos o gráfico:

Com base no gráfico, as previsões parecem acompanhar razoavelmente bem a sazonalidade dos dados originais, embora fique abaixo ou acima dos valores reais em alguns períodos. A maior discrepância é percebida para o mês de abril, onde o valor predito foi de aprox. R$ 126.175,10 , contra R$ 144.770,00 de valor verdadeiramente vendido. Apesar disso, para os demais meses a diferença é bem menos abrupta.
Assim sendo, as previsões parecem ser adequadas.

Item a)

Faça um gráfico da série.

Item b)

Observando o gráfico, tente adivinhar um valor aproximado para a autocorrelação em h = 1.

Olhando para o gráfico, a autocorrelação em h=1 parece ser forte e negativa, como algo em torno de -0.7.
Isso se dá porque quanto maior (menor) o valor no instante t, menor (maior) é o valor em t+1, quase como se fossem “saltos simétricos”.

Item c)

Faça um gráfico de \(x_t\) e \(x_{t+1}\) e novamente tente adivinhar um valor para ρ(1).

De fato, parece haver a estrutura de uma reta com inclinação negativa. No entanto, esse gráfico já dá a ideia de uma correlação moderada, e não forte como “adivinhado” anteriormente. Com base no gráfico atual, um possível valor para a autocorrelação em h=1 seria algo em torno de -0.5.

Item d)

Calcule \(\hat{\rho}(1)\).

O estimador da autocorrelação será dado por:

\[ \hat{\rho}(1) = \frac{\hat{\gamma}(1)}{\hat{\gamma}(0)} = \frac{\sum_{t=1}^{n-1} (X_t - \bar{X})(X_{t+1} - \bar{X})}{\sum_{t=1}^{n}(X_t - \bar{X})^2} \]

xbarra=mean(valor)
n=length(valor)

t=1
numerador=0
denominador=0

while(t<=I(n-1)){
  numerador=numerador+(valor[t]-xbarra)*(valor[t+1]-xbarra)
  denominador=denominador+(valor[t]-xbarra)^2
  t=t+1
}

denominador=denominador+(valor[t]-xbarra)^2

rho_chapeu=numerador/denominador

O valor calculado para \(\hat{\rho}(1)\) foi -0.5487805.

Série Ozônio

result_oz=data.frame(h=0:4,gamma=NA,rho=NA)

result_oz[result_oz$h==0,"gamma"]=autocov(ozonio$Ozonio,0)
result_oz[result_oz$h==1,"gamma"]=autocov(ozonio$Ozonio,1)
result_oz[result_oz$h==2,"gamma"]=autocov(ozonio$Ozonio,2)
result_oz[result_oz$h==3,"gamma"]=autocov(ozonio$Ozonio,3)
result_oz[result_oz$h==4,"gamma"]=autocov(ozonio$Ozonio,4)

result_oz[result_oz$h==0,"rho"]=autocorr(ozonio$Ozonio,0)
result_oz[result_oz$h==1,"rho"]=autocorr(ozonio$Ozonio,1)
result_oz[result_oz$h==2,"rho"]=autocorr(ozonio$Ozonio,2)
result_oz[result_oz$h==3,"rho"]=autocorr(ozonio$Ozonio,3)
result_oz[result_oz$h==4,"rho"]=autocorr(ozonio$Ozonio,4)

rownames(result_oz)=0:4

Série Energia

result_en=data.frame(h=0:4,gamma=NA,rho=NA)

result_en[result_en$h==0,"gamma"]=autocov(energia$Energia,0)
result_en[result_en$h==1,"gamma"]=autocov(energia$Energia,1)
result_en[result_en$h==2,"gamma"]=autocov(energia$Energia,2)
result_en[result_en$h==3,"gamma"]=autocov(energia$Energia,3)
result_en[result_en$h==4,"gamma"]=autocov(energia$Energia,4)

result_en[result_en$h==0,"rho"]=autocorr(energia$Energia,0)
result_en[result_en$h==1,"rho"]=autocorr(energia$Energia,1)
result_en[result_en$h==2,"rho"]=autocorr(energia$Energia,2)
result_en[result_en$h==3,"rho"]=autocorr(energia$Energia,3)
result_en[result_en$h==4,"rho"]=autocorr(energia$Energia,4)

rownames(result_en)=0:4
result_en

  h      gamma       rho
0 0 2895122894 1.0000000
1 1 2754466348 0.9514160
2 2 2656106971 0.9174419
3 3 2559293421 0.8840017
4 4 2456155749 0.8483770

Resultados de autocovariância e autocorrelação para a série Ozônio:

result_oz

  h      gamma        rho
0 0  4.3114108  1.0000000
1 1  2.9744369  0.6898987
2 2  1.5887360  0.3684956
3 3 -0.1794922 -0.0416319
4 4 -1.7265890 -0.4004696

Resultados de autocovariância e autocorrelação para a série Energia:

result_en

  h      gamma       rho
0 0 2895122894 1.0000000
1 1 2754466348 0.9514160
2 2 2656106971 0.9174419
3 3 2559293421 0.8840017
4 4 2456155749 0.8483770

Gráficos

Vejamos os gráficos para a série Ozônio:

Para a série Ozônio, a presença de sazonalidade fica claramente evidenciada no gráfico, com a série oscilando de forma bem comportada em torno de sua média (5.0794), sem a presença de tendência. Isso nos sugere que a quantidade de Ozônio varia de acordo com os meses do ano. A sazonalidade também é capturada pelo gráfico da função de autocorrelação, que mostra alternância entre correlações positivas e negativas de forma cíclica e bem comportada.

Agora vamos aos gráficos para a série Energia:

Já com respeito à série de Energia, verificamos uma estrutura de tendência crescente ao longo dos anos/meses, sem indicação da existência de sazonalidade. Isso também é confirmado pelo gráfico da função de correlação, que decresce à medida que h aumenta, porém sem apresentar oscilações.