Integrantes do Grupo:

Introdução e Tratamento dos Dados

Usaremos o conjunto de pacotes "tidyverse" para realizar algumas operações, sobretudo utilizando recursos gráficos.

# ---- Importação das Bibliotecas ----
library("tidyverse")
## Warning: package 'ggplot2' was built under R version 4.0.5
## Warning: package 'readr' was built under R version 4.0.4
library("tidyquant")
library("Quandl")
library("timetk")
## Warning: package 'timetk' was built under R version 4.0.5
library("devtools")
library("quantmod")
library("readxl")
library("GetQuandlData")
library("AER")
library("strucchange")
# ---- Importação dos Dados ----

# Cotações
ativos <- c("VALE3.SA", "^BVSP")
inicio <- as.Date("2010-01-04")
fim <- as.Date("2020-12-30")

cotacao <- ativos %>% 
    tq_get(get = "stock.prices", from = inicio, to = fim) %>% 
    mutate(adjusted = na.locf(adjusted))

# Retornos

df_retornos <- cotacao %>% 
    group_by(symbol) %>% 
    select(date, adjusted) %>% 
    tq_transmute(select = adjusted,
                 mutate_fun = periodReturn,
                 period = "daily",
                 type = "arithmetic",
                 col_rename = "retorno") %>% 
    spread(key = symbol, value = retorno)

# ---- Importando a taxa SELIC da API do BACEN (de 2010 até hoje) ----

Quandl.api_key('TC1ow5j6G7s4SFHTzgDz') # Acrescentando a chave da API - Comando necessário pra acessar o Quandl
selic_diaria <-  as.data.frame(Quandl("BCB/11", type = "xts")/100) # Importando a serie do selic do Bacen
df_selic <- selic_diaria %>% 
    rownames_to_column() %>% 
    tibble() %>% 
    rename("date" = rowname, "selic" = V1) %>% 
    mutate(selic = na.locf(selic))

df_selic$date <- as.Date(df_selic$date)

# acréscimo na base de dados
df_retornos <- df_retornos %>% 
    left_join(df_selic, by = "date")


Respostas


# ---- Montando regressão do CAPM ---- 
df_capm <- df_retornos %>% 
    transmute("date" = date, "Rp-Rf" = VALE3.SA - selic, "Rm-Rf" = `^BVSP` - selic)

capm1 <- lm(df_capm$`Rp-Rf` ~ df_capm$`Rm-Rf`)
summary(capm1)
## 
## Call:
## lm(formula = df_capm$`Rp-Rf` ~ df_capm$`Rm-Rf`)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.221786 -0.010506 -0.000477  0.010061  0.114318 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     0.0003861  0.0003983   0.969    0.332    
## df_capm$`Rm-Rf` 1.0576230  0.0250455  42.228   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02078 on 2719 degrees of freedom
##   (6 observations deleted due to missingness)
## Multiple R-squared:  0.3961, Adjusted R-squared:  0.3959 
## F-statistic:  1783 on 1 and 2719 DF,  p-value: < 2.2e-16
# ---- Montando regressão de Fama-French ----

# (adaptar os arquivos diretamente no Excel, manipulando a data para o formato YYYY-MM-DD)
# Posso usar a função as.Date(paste(year, month, day, sep = "-")) também
MKT_Factor <- read_excel("D:/Dropbox/UFMG/Programação/R/Métodos/Ativ Prát 3/Market_Factor.xls")
## New names:
## * `` -> ...6
HML_Factor <- read_excel("D:/Dropbox/UFMG/Programação/R/Métodos/Ativ Prát 3/HML_Factor.xls")
## New names:
## * `` -> ...6
SMB_Factor <- read_excel("D:/Dropbox/UFMG/Programação/R/Métodos/Ativ Prát 3/SMB_Factor.xls")
## New names:
## * `` -> ...6
HML <- HML_Factor %>% select(date, HML)
HML$date <- as.Date(HML$date)

MKT <- MKT_Factor %>% select(date, Rm_minus_Rf) %>% rename("MKT" = Rm_minus_Rf)
MKT$date <- as.Date(MKT$date)

SMB <- SMB_Factor %>% select(date, SMB)
SMB$date <- as.Date(SMB$date)

QUESTÃO 1:

Montando gráfico com os retornos da ação escolhida

# Gráfico CAPM

graf_1 <- df_capm %>% 
    ggplot(aes(x = `Rm-Rf`, y = `Rp-Rf`)) +
    geom_point() +
    labs(x = "Prêmio de Risco do Mercado",
         y = "Prêmio de Risco do Ativo",
         title = "CAPM",
         subtitle = "Empresa: VALE S.A",
         caption = "Dados coletados do Yahoo Finance e da API do Quandl") +
    geom_hline(yintercept = 0, linetype = "dashed") +
    geom_vline(xintercept = 0, linetype = "dashed") +
    geom_smooth(method = "lm", color = "red") + 
    theme_bw()

graf_1
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 6 rows containing non-finite values (stat_smooth).
## Warning: Removed 6 rows containing missing values (geom_point).

graf_2 <- MKT %>% 
    ggplot(aes(x = date, y = MKT)) +
    geom_point() +
    labs(x = "",
         y = "MKT",
         title = "Fator MKT",
         caption = "Dados coletados do site NEFIN-USP") +
    geom_hline(yintercept = 0, linetype = "dashed") +
    theme_bw()

graf_2

graf_3 <- HML %>% 
    ggplot(aes(x = date, y = HML)) +
    geom_point() +
    labs(x = "",
         y = "HML",
         title = "Fator HML",
         caption = "Dados coletados do site NEFIN-USP") +
    geom_hline(yintercept = 0, linetype = "dashed") +
    theme_bw()

graf_3

graf_4<- SMB %>% 
    ggplot(aes(x = date, y = SMB)) +
    geom_point() +
    labs(x = "",
         y = "SMB",
         title = "Fator SMB",
         caption = "Dados coletados do site NEFIN-USP") +
    geom_hline(yintercept = 0, linetype = "dashed") +
    theme_bw()

graf_4

SMB <- SMB_Factor %>% select(date, SMB)
SMB$date <- as.Date(SMB$date)

# Gráfico Retornos 

graf_5 <- df_retornos %>% 
    ggplot(aes(x = date, y = VALE3.SA)) +
    geom_line() +
    geom_hline(yintercept = 0, linetype = "dashed") +
    scale_y_continuous(labels = scales::percent) +
    labs(title = "Retornos Mensais",
         subtitle = "Empresa: Petrobrás S.A",
         y = "Retornos Mensais",
         x = "") + 
    theme_tq() + 
    scale_fill_tq()

graf_5

# ---- Montando regressão de três fatores de Fama-French ----
fatores <- df_capm %>% 
    left_join(HML, by = "date") %>% 
    left_join(SMB, by = "date") %>% 
    left_join(MKT, by = "date") %>% 
    select(-`Rm-Rf`)
    
ffrench1 <- lm(`Rp-Rf` ~ MKT + HML + SMB, data = fatores)
summary(ffrench1)
## 
## Call:
## lm(formula = `Rp-Rf` ~ MKT + HML + SMB, data = fatores)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.230730 -0.010448 -0.000378  0.009950  0.133562 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.0002042  0.0003985   0.512    0.608    
## MKT          0.9891863  0.0294807  33.554  < 2e-16 ***
## HML          0.7748240  0.0551971  14.037  < 2e-16 ***
## SMB         -0.4052043  0.0507856  -7.979 2.16e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02077 on 2714 degrees of freedom
##   (9 observations deleted due to missingness)
## Multiple R-squared:  0.3976, Adjusted R-squared:  0.3969 
## F-statistic: 597.1 on 3 and 2714 DF,  p-value: < 2.2e-16

QUESTÃO 2:

O teste F verifica se o modelo restrito é valido ou não. Quando a estatística F for alta, assim como obtivemos no resultado do teste, terá um p-valor significativo. Isso demonstra que a média condicional (que é estimada pelo modelo) é superior a média incondicional, logo, rejeita-se a hipótese numa que os dois modelos tem a mesma eficácia.

# Teste F entre os dois modelos
anova(capm1,ffrench1)
## Warning in anova.lmlist(object, ...): models with response '"Rp-Rf"' removed
## because response differs from model 1
## Analysis of Variance Table
## 
## Response: df_capm$`Rp-Rf`
##                   Df  Sum Sq Mean Sq F value    Pr(>F)    
## df_capm$`Rm-Rf`    1 0.76968 0.76968  1783.2 < 2.2e-16 ***
## Residuals       2719 1.17359 0.00043                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

QUESTÃO 3:

Valores médios para VIF superiores a 5 s?o um indício de multicolinearidade, sendo assim não há multicolineariade nos testes pois todos os valores são menores que 5.

# VIF do modelo Fama-French
vif(ffrench1)
##      MKT      HML      SMB 
## 1.101118 1.221621 1.141517

QUESTÃO 4:

Foi obtido estatísticas de teste altas, de forma que o p-valor foi significativo para os dois modelos. Isso indica que há uma forma funcional ou especificação alternativa que melhor valide o modelo. Fazendo o teste para os regressores também identificamos p-valores significativos, onde devemos verificar se as variáveis independentes são melhores ajustadas ao modelo se tratadas de forma polinomial.

# Teste RESET em cada modelo 

resettest(capm1, power = c(2,3,4))
## 
##  RESET test
## 
## data:  capm1
## RESET = 6.5551, df1 = 3, df2 = 2716, p-value = 0.0002056
resettest(ffrench1, power = c(2,3,4))
## 
##  RESET test
## 
## data:  ffrench1
## RESET = 7.153, df1 = 3, df2 = 2711, p-value = 8.784e-05
# Avalia o teste como potências dos regressores, e não potências de todas as variáveis em conjunto

resettest(ffrench1, power = c(2,3,4), type = "regressor") 
## 
##  RESET test
## 
## data:  ffrench1
## RESET = 6.7626, df1 = 9, df2 = 2705, p-value = 1.191e-09

QUESTÃO 5:

# Teste de Chow (quebra estrutural)

sctest(df_capm$`Rp-Rf` ~ df_capm$`Rm-Rf`, data = df_capm, type = "Chow")
## 
##  Chow test
## 
## data:  df_capm$`Rp-Rf` ~ df_capm$`Rm-Rf`
## F = 3.2026, p-value = 0.04081
sctest(`Rp-Rf` ~ MKT + HML + SMB, data = fatores, type = "Chow")
## 
##  Chow test
## 
## data:  `Rp-Rf` ~ MKT + HML + SMB
## F = 6.64, p-value = 2.579e-05
# Gráfico do teste F para subgrupos da amostra
capm_f <- Fstats(df_capm$`Rp-Rf` ~ df_capm$`Rm-Rf`, data = df_capm, from = 0.1, to = 0.9)
plot(capm_f, main = "CAPM Estatísticas F") # Gráfico do F

plot(capm_f, pval = T, main = "CAPM p-valores") # Gráfico do p-valor

ffrench1_f <- Fstats(`Rp-Rf` ~ MKT + HML + SMB, data = fatores, from = 0.1, to = 0.9)
plot(ffrench1_f, main = "Fama-French Estatísticas F") # Gráfico do F

plot(ffrench1_f, pval = T, main = "Fama-French p-valores") # Gráfico do p-valor

QUESTÃO 6:

O modelo de Fama-French é preferível, dado o teste F;

O teste RESET indica que há formas funcionais dos REGRESSORES que se ajustam melhor ao modelo, tanto no CAPM quanto no Fama-French. Isso indica que o fator “prêmio de risco do mercado” pode ser ajustado para a obtenção de melhores resultados;

O Teste de Chow indica que há quebras estruturais em subgrupos das amostras, ou seja, outra regressão pode ser aplicada em algum ponto do modelo para que este seja melhor ajustado.