Primeiro, definimos nossa base de dados e estabelecemos o modelo capm1 e ffrench.

dados <- read_excel("Atividade 3.xlsx")
attach(dados)
## # A tibble: 6 x 9
##   Data     Market      HML      SMB      WML Risk_free    VALE3   Size1      BM1
##   <chr>     <dbl>    <dbl>    <dbl>    <dbl>     <dbl>    <dbl>   <dbl>    <dbl>
## 1 04/01~  0.0214   3.45e-3 -0.00544 -0.00640  0.000329  0.0401  1.49e-2  1.38e-2
## 2 05/01~  0.00133  4.97e-3 -0.00411 -0.00674  0.000330  0.00925 5.38e-5 -3.38e-3
## 3 06/01~  0.00587  4.73e-3  0.0144   0.00624  0.000329  0.0211  1.87e-2  1.05e-2
## 4 07/01~ -0.00295  5.33e-3  0.00984  0.00367  0.000329  0.00419 7.09e-3  8.72e-4
## 5 08/01~ -0.00184 -2.92e-5  0.00961 -0.00294  0.000329  0.00984 8.57e-3  4.46e-3
## 6 11/01~  0.00147  1.23e-2  0.0100   0.00545  0.000329 -0.00295 1.04e-2  9.11e-4
vale_rf <- VALE3 - Risk_free
capm1 <- lm(vale_rf ~ Market)
## 
## Call:
## lm(formula = vale_rf ~ Market)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.228209 -0.010489 -0.000388  0.010256  0.119028 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.0003088  0.0004082   0.757    0.449    
## Market      1.1296982  0.0288118  39.210   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.0213 on 2721 degrees of freedom
## Multiple R-squared:  0.361,  Adjusted R-squared:  0.3608 
## F-statistic:  1537 on 1 and 2721 DF,  p-value: < 2.2e-16
ffrench <- update(capm1, formula=.~.+ SMB + HML)
## 
## Call:
## lm(formula = vale_rf ~ Market + SMB + HML)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.230339 -0.010460 -0.000362  0.009916  0.133130 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.0002566  0.0003928   0.653    0.514    
## Market       0.9987656  0.0290854  34.339  < 2e-16 ***
## SMB         -0.4128445  0.0500822  -8.243 2.57e-16 ***
## HML          0.7874732  0.0544349  14.466  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02049 on 2719 degrees of freedom
## Multiple R-squared:  0.4091, Adjusted R-squared:  0.4084 
## F-statistic: 627.4 on 3 and 2719 DF,  p-value: < 2.2e-16

Questão 1

Teste Jarque-Bera

jarqueberaTest(resid(capm1))
## 
## Title:
##  Jarque - Bera Normalality Test
## 
## Test Results:
##   STATISTIC:
##     X-squared: 7379.8273
##   P VALUE:
##     Asymptotic p Value: < 2.2e-16 
## 
## Description:
##  Thu Feb 11 18:22:02 2021 by user: Pedro
jarqueberaTest(resid(ffrench))
## 
## Title:
##  Jarque - Bera Normalality Test
## 
## Test Results:
##   STATISTIC:
##     X-squared: 10340.9115
##   P VALUE:
##     Asymptotic p Value: < 2.2e-16 
## 
## Description:
##  Thu Feb 11 18:22:02 2021 by user: Pedro

Como pode ser observado, a hipótese nula de normalidade dos resíduos pode ser rejeitada, dado que o valor da estatística qui-quadrado foi bastante expressivo e, consequentemente, o p-valor foi bastante pequeno.

Observações Influentes

par(mfrow=c(2,3))
plot(capm1, which=1:6)

summary(influence.measures(capm1))
par(mfrow=c(2,3))
plot(ffrench, which=1:6)

summary(influence.measures(ffrench))

Devido ao grande número de observações que atendiam a algum critério na definição de dado discrepante (centenas), foi necessário suprimir o output dessa função para que ela não ocupasse um espaço muito grande no documento.

Entretanto, como fica claro pelos gráficos e pelo tamanho da resposta para a função “summary(influence.measures())”, existem diversos outliers em nossa amostra.

Heterocedasticidade

bptest(capm1)
## 
##  studentized Breusch-Pagan test
## 
## data:  capm1
## BP = 0.17028, df = 1, p-value = 0.6799
bptest(ffrench)
## 
##  studentized Breusch-Pagan test
## 
## data:  ffrench
## BP = 12.246, df = 3, p-value = 0.006585

Pelo teste Breusch-Pagan não podemos rejeitar a hipótese nula de homocedasticidade do modelo capm1. Entretanto, para o modelo ffrench, podemos rejeitar a hipótese nula de homocedasticidade. O fato do modelo restrito (capm1) ser homocedástico enquanto o modelo irrestrito (ffrench) é heterocedástico, nos indica que alguma variável aficionada ao modelo ffrench (SMB, HML) pode ser a causadora da heterocedasticidade.

Correlação Serial

dwtest(capm1)
## 
##  Durbin-Watson test
## 
## data:  capm1
## DW = 1.8782, p-value = 0.000744
## alternative hypothesis: true autocorrelation is greater than 0
dwtest(ffrench)
## 
##  Durbin-Watson test
## 
## data:  ffrench
## DW = 1.9002, p-value = 0.004581
## alternative hypothesis: true autocorrelation is greater than 0
bgtest(capm1)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  capm1
## LM test = 10.048, df = 1, p-value = 0.001525
bgtest(ffrench)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  ffrench
## LM test = 6.7677, df = 1, p-value = 0.009282

Como o resultado do teste de Durbin-Watson retornou valores bastante próximos de 2 em ambos os modelos, podemos afirmar que não existe correlação serial de primeira ordem. Entretanto, o teste de Breusch-Godfrey nos indica que existe sim correlação serial em ordens superiores.

Questão 2

Por “corrigir para a presença de eventuais outliers” entendemos suprimir da amostra elementos que podem ser considerados como discrepantes.

detach(dados)
out_capm1 <- which(apply(influence.measures(capm1)$is.inf, 1, any))
out_ffrench <- which(apply(influence.measures(ffrench)$is.inf, 1, any))

Assim, nossas novas regressões, sem outliers, são:

capm1_out <- lm(I(VALE3 - Risk_free) ~ Market, data = dados[-out_capm1,])
## 
## Call:
## lm(formula = I(VALE3 - Risk_free) ~ Market, data = dados[-out_capm1, 
##     ])
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.047049 -0.009599  0.000006  0.009831  0.046910 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0002421  0.0003238  -0.748    0.455    
## Market       1.1423718  0.0305823  37.354   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.0162 on 2503 degrees of freedom
## Multiple R-squared:  0.3579, Adjusted R-squared:  0.3577 
## F-statistic:  1395 on 1 and 2503 DF,  p-value: < 2.2e-16
ffrench_out <- lm(I(VALE3 - Risk_free) ~ Market + HML + SMB, data = dados[-out_ffrench,])
## 
## Call:
## lm(formula = I(VALE3 - Risk_free) ~ Market + HML + SMB, data = dados[-out_ffrench, 
##     ])
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.046559 -0.009501  0.000095  0.009496  0.045707 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0002448  0.0003091  -0.792    0.429    
## Market       0.9067059  0.0308054  29.433  < 2e-16 ***
## HML          0.7611039  0.0512794  14.842  < 2e-16 ***
## SMB         -0.2874659  0.0474529  -6.058 1.59e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.0154 on 2477 degrees of freedom
## Multiple R-squared:  0.416,  Adjusted R-squared:  0.4153 
## F-statistic: 588.1 on 3 and 2477 DF,  p-value: < 2.2e-16

Teste Jarque-Bera

jarqueberaTest(resid(capm1_out))
## 
## Title:
##  Jarque - Bera Normalality Test
## 
## Test Results:
##   STATISTIC:
##     X-squared: 2.5437
##   P VALUE:
##     Asymptotic p Value: 0.2803 
## 
## Description:
##  Thu Feb 11 18:22:08 2021 by user: Pedro
jarqueberaTest(resid(ffrench_out))
## 
## Title:
##  Jarque - Bera Normalality Test
## 
## Test Results:
##   STATISTIC:
##     X-squared: 1.8747
##   P VALUE:
##     Asymptotic p Value: 0.3917 
## 
## Description:
##  Thu Feb 11 18:22:08 2021 by user: Pedro

Como podemos observar, após a retirda dos outliers, temos que os resíduos são normais, diferentemente do modelo com os outliers.

Heterocedasticidade

bptest(capm1_out)
## 
##  studentized Breusch-Pagan test
## 
## data:  capm1_out
## BP = 7.6772, df = 1, p-value = 0.005592
bptest(ffrench_out)
## 
##  studentized Breusch-Pagan test
## 
## data:  ffrench_out
## BP = 3.5579, df = 3, p-value = 0.3133

Como é possível observar, houve uma inversão. O modelo CAPM com outlier é homocedástico, enquanto o modelo sem os outliers é heterocedástico. O modelo Fama & French com outliers é heterocedástico, enquanto sem outliers passa a ser homocedástico.

Correlação Serial

dwtest(capm1_out)
## 
##  Durbin-Watson test
## 
## data:  capm1_out
## DW = 1.9833, p-value = 0.338
## alternative hypothesis: true autocorrelation is greater than 0
dwtest(ffrench_out)
## 
##  Durbin-Watson test
## 
## data:  ffrench_out
## DW = 1.9416, p-value = 0.07251
## alternative hypothesis: true autocorrelation is greater than 0
bgtest(capm1_out)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  capm1_out
## LM test = 0.16238, df = 1, p-value = 0.687
bgtest(ffrench_out)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  ffrench_out
## LM test = 2.0778, df = 1, p-value = 0.1495

Através desses testes podemos observar que, após a retirada dos outliers, eliminamos a correlação serial.

Questão 3

Podemos observar melhoras significativas nas estatísticas de teste retirando os outliers, entretanto, esse procedimento não é muito aconselhável, dado que esses elementos retirados são parte do PGD. Segundo Draper & Smith 1:

“A rejeição automática de dados discrepantes nem sempre é um procedimento sensato. Às vezes o dado discrepante está fornecendo informações que outros não podem fornecer, devido ao fato de que ele decorre de uma combinação incomum de circunstância que podem ser de interesse vital e exigem mais investigação, e não rejeição. Como regra geral, os dados discrepantes só devem ser rejeitados se puderem ser atribuídos a causas como erros de registro das observações ou problemas de precisão de aparelhos. Caso contrário, merecem investigação cuidadosa.”

Ou seja, não é aconselhável retirar os outliers da amostra. Especialmente, no caso da VALE3, estaríamos eliminando justamente as amostras que podem ser mais significativas para o investidor - uma queda (alta) incomum de 10% em um dia é mais relevante para o investidor do que pequenas variações de 0,1%. Os efeitos dessa retirada podem ser vistos facilmente quando analisamos o retorno anualizado e desvio padrão anualizado com e sem os outliers para o CAPM:

((prod(1 + dados[-out_capm1,]$VALE3))^(252/length(dados[-out_capm1,]$VALE3))-1)*100
## [1] 4.283976
sd(dados[-out_capm1,]$VALE3)*sqrt(252)*100
## [1] 32.07401
((prod(1 + dados$VALE3))^(252/length(dados$VALE3))-1)*100
## [1] 10.64591
sd(dados$VALE3)*sqrt(252)*100
## [1] 42.28796

Como é possível observar, houve uma mudança significativa tanto no padrão de comportamento do retorno acumulado, quanto no risco da amostra. Fazer inferências com essa amostra reduzida pode ser bastante problemático, dado que houveram alterações substanciais no padrão de comportamento dos retornos após a exclusão dos outliers. É interessante observar também que a retirada dos outliers reduziu o retorno total do ativo, indicando que boa parte dos outliers eram dias em que houve um retorno positivo anormal. Isso pode ser um indicativo de que temos assimetria positiva na amostra.

Correção de problemas

Não normalidade dos resíduos, outliers e heterocedasticidade.

No modelo original percebemos que temos múltiplas amostras influentes (outliers), não normalidade dos resíduos, correlação serial e heterocedasticidade. Visando solucionar esses problemas iremos realizar uma regressão robusta (Least Absolute Deviations) e adicionar defassagens dos retornos para solucionar (ou mitigar) o probelma de correlação serial. A vantagem dessa regressão robusta é que ela não faz nenhuma suposição acerca da normalidade dos resíduos e dá pouca importância para outliers, além de ser robusto para a heterocedasticidade 2.

capm2 <- rq(I(VALE3 - Risk_free) ~ Market, data = dados)
## 
## Call: rq(formula = I(VALE3 - Risk_free) ~ Market, data = dados)
## 
## tau: [1] 0.5
## 
## Coefficients:
##             Value    Std. Error t value  Pr(>|t|)
## (Intercept) -0.00010  0.00037   -0.27805  0.78099
## Market       1.10670  0.02621   42.21979  0.00000
ffrench2 <- rq(vale_rf ~ Market + SMB + HML, data = dados)
summary(ffrench2)
## 
## Call: rq(formula = vale_rf ~ Market + SMB + HML, data = dados)
## 
## tau: [1] 0.5
## 
## Coefficients:
##             Value    Std. Error t value  Pr(>|t|)
## (Intercept) -0.00008  0.00034   -0.24331  0.80779
## Market       0.93871  0.01487   63.14580  0.00000
## SMB         -0.34252  0.03939   -8.69458  0.00000
## HML          0.70586  0.03938   17.92475  0.00000

Questão 4

Para fins de inferência a regressão robusta (LAD) deve ser utilizada. Com essa regressão solucionamos os problemas relacionados a outliers na amostra, a heterocedasticidade e não normalidade dos resíduos. Segundo Hou, et al 3, a aplicação do teste Breusch-Godfrey “pode resultar em grandes distorções, especialmente para o primeiro ou último quartil” quando o utilizamos para regressões quantílicas em time series. O autor propõe um teste alternativo (QF Test), entretanto esse teste foge dos conhecimentos dos autores do presente trabalho. Assim, assumimos a limitação do modelo, mas continuamos defendendo sua superioridade frente ao tradicional MQO, nesse caso, pelas razões já expostas.

acf(na.omit(capm2$residuals), plot = T)

Além disso, a correlação serial, nesse caso, é uma preocupação menor, dado que, como pode ser visto pelo gráfico, o intervalo de confiança é desrespeitado em poucos casos e por uma margem pequena.

  1. Gujarati, p. 496↩

  2. James L Powell, Least absolute deviations estimation for the censored regression model, Journal of Econometrics, Volume 25, Issue 3, 1984, Pages 303-325, ISSN 0304-4076, https://doi.org/10.1016/0304-4076(84)90004-6.↩

  3. Testing for Autocorrleation in Quantile Regression Models; Lijuan Huo, Tae-Hwan Kim, Yunmi Kim and Dong Jin Lee; 2014↩