Pacotes Utilizados, Download e Tratamento dos Dados

Foram utilizados os pacotes rugarch, ggplot2, quantmod, BatchGetSymbols, kableExtra, lubridate, gridExtra, tseries, reshape2 e dplyr.

O código abaixo realiza o download dos dados da série diária do Ibovespa a partir da data definida em “start” até a data definida em “end”, através da função BatchGetSymbols.

start<-as.Date("2000-01-01")
end<-as.Date("2020-05-29")
ticker<-c("^BVSP")

bvsp = BatchGetSymbols(ticker, first.date = start, last.date = end)
## 
## Running BatchGetSymbols for:
##    tickers =^BVSP
##    Downloading data for benchmark ticker
## ^GSPC | yahoo (1|1) | Found cache file
## ^BVSP | yahoo (1|1) | Found cache file - Got 97% of valid prices | Got it!
bvsp<-bvsp$df.tickers
bvsp$ref.date<-as.Date(bvsp$ref.date)

bvsp_xts<-na.omit(xts(bvsp[,4],order.by = bvsp$ref.date))
dates<-rownames(bvsp_xts)
colnames(bvsp_xts)<-"IBOV"

Abaixo é apresentado um gráfico da série obtida, com 5048 observações diárias, iniciada em 2000-01-01 e com término em 2020-05-29.

Análise dos Retornos Diários do Ibovespa

Obtemos, então, os retornos diários da série do Ibovespa e elaboramos o gráfico 2, abaixo, que apresenta o comportamento dos retornos diários do Ibovespa. Nele, é possível perceber a presença de aglomerados de volatilidade e média aparentemente constante em zero, ilustrada pela linha vermelha.

Outro fato estilizado que pode ser verificado através dos gráficos 3 e 4, abaixo, é a distribuição dos retornos não ser Normal-Gaussiana: as séries dos retornos diários do Ibovespa apresentam assimetria negativa: retornos negativos de maior valor ocorrem mais frequentemente que retornos positivos mais altos. Também apresentam caudas mais grossas, indicando maior probabilidade de eventos extremos do que o previsto por uma distribuição Normal, conforme fica deveras evidente no segundo gráfico. Nestes dois últimos gráficos, a área rosa corresponde a uma distribuição Normal com média e desvio padrão extraídos da série dos retornos diários do Ibovespa, enquanto a área azul corresponde à distribuição dos retornos do Ibovespa.

## Don't know how to automatically pick scale for object of type xts/zoo. Defaulting to continuous.

## Don't know how to automatically pick scale for object of type xts/zoo. Defaulting to continuous.

Testamos a hipótese de não normalidade da série dos retornos através do teste Jarque-Bera:

## 
##  Jarque Bera Test
## 
## data:  bvsp_ret
## X-squared = 8241, df = 2, p-value < 2.2e-16

Uma vez que \(p-valor<0.05\), temos evidências suficientes para rejeitar a hipótese nula de que os retornos são normalmente distribuídos. Assumir que os retornos seguem uma distribuição normal pode levar a uma subestimação do risco. Portanto, assumindo que a distribuição dos retornos é mais próxima de uma t-student, simétrica, porém com caudas mais pesadas, pode ser uma premissa adequada.

Estimação do GARCH(1,1)

Uma vez que os retornos apresentam heterocedasticidade condicional, evidenciada pelos aglomerados de volatilidade, e distribuição com caudas pesadas, vamos estimar abaixo um modelo GARCH(1,1):

garch.model<-ugarchspec(variance.model = list(model = 'sGARCH' , garchOrder = c(1 , 1)) , 
                        mean.model = list(armaOrder = c(0,0),include.mean=FALSE),
                        distribution.model = "std")
fit.garch <- ugarchfit(spec = garch.model , data = bvsp_ret)

Os coeficientes encontrados são apresentados abaixo. Cabe ressaltar que todos os parâmetros são estatisticamente significantes.

##            Estimate   Std. Error   t value     Pr(>|t|)
## omega  5.137252e-06 1.985375e-06  2.587548 9.666180e-03
## alpha1 6.837927e-02 7.325239e-03  9.334750 0.000000e+00
## beta1  9.147711e-01 9.251920e-03 98.873653 0.000000e+00
## shape  1.036783e+01 1.319231e+00  7.858992 3.774758e-15

O modelo GARCH(1,1) estimado é: \[\sigma_t^2=0.00001+0.06838r_{t-1}^2+0.91477\sigma_{t-1}^2\] Para estimar a variância incondicional, \(V_L\), calculamos: \[V_L=\frac{\omega}{1-(\alpha+\beta)}\]

omega<-fit.garch@fit$coef[1]
alpha<-fit.garch@fit$coef[2]
beta<-fit.garch@fit$coef[3]
Vl<-round((omega/(1-alpha-beta))*100,2)
Vla<-Vl*252

0u seja, a variância de longo prazo é de 0.03% ao dia, ou 7.56% ao ano. Para efeito de comparação, a variância da série dos retornos do Ibovespa no período analisado é de 0.03% ao dia, ou 7.56% ao ano.Abaixo, podemos constatar a variância condicional do modelo GARCH(1,1) com a série dos retornos do Ibovespa ao quadrado e a variância incondicional da série:

De forma semelhante, podemos constatar abaixo a volatilidade condicional do modelo GARCH(1,1) com a série dos retornos do Ibovespa e a volatilidade incondicional da série:

Previsão da Volatilidade Através do GARCH(1,1)

Com o intuito de previsão, abaixo são apresentados os valores ao quadrado da série dos retornos do Ibovespa e sua variância condicional nos últimos 30 dias e a previsão da para os próximos 30 dias da variância condicional do modelo GARCH(1,1) estimado.

h<-30
b<-30
ugfore<-ugarchforecast(fit.garch,n.ahead=h)
ug_f<-ugfore@forecast$sigmaFor

fit_var_tail<-c(tail(fit_var,b),rep(NA,h))
fit_res2_tail<-c(tail(fit_res2,b),rep(NA,h))
Vl_f<-c(rep(Vl/100,h+b))
Vl_f<-as.numeric(Vl_f)
forecast<-c(rep(NA,b),(ug_f)^2)
index<-(-b:h)
index<-index[-((b+h)/2+1)]