Introdução

Carregando Pacotes e Dados


Carregando Pacotes

# Pacotes para manipulação
library(tidyverse)
library(scales)
library(tidyquant)
library(lubridate)
library(writexl)

# Pacotes para modelagem
library(PerformanceAnalytics)
library(xts)
library(forecast)

# Pacotes para dados
library(rbcb)
library(rb3)
library(ggthemes)
library(rvest) # webscraping
library(zoo)


Carregando os Dados

# nomes_ativos <- "https://br.tradingview.com/markets/stocks-brazil/market-movers-large-cap/" %>% 
#   read_html() %>%   
#   html_elements(".tickerName-GrtoTeat") %>% 
#   html_text() %>%
#   str_c(".SA")
# 
# ativos = tq_get(x = nomes_ativos,
#                 from = '2013-01-01',
#                 to = '2023-12-31')
# 
# ibov = tq_get(x = '^BVSP',
#               from = '2013-01-01',
#               to = '2023-12-31')
# 
# ativos_filtrado = ativos %>% 
#   dplyr::mutate(month = month(date),
#                 year = year(date),
#                 day = day(date)) %>% 
#   dplyr::group_by(symbol) %>% 
#   dplyr::filter(all(year >= 2013)) %>%
#   dplyr::ungroup() %>% 
#   dplyr::select(-c(month, year, day)) %>% 
#   dplyr::mutate(date = make_date(year(date),
#                                  month(date),
#                                  day(date))) %>% 
#   dplyr::arrange(date) %>% 
#   dplyr::rename(acao = adjusted, ticker = symbol) %>% 
#   dplyr::select(ticker, date, acao)
# 
# ativos_modelagem = ativos_filtrado %>% 
#   tidyr::pivot_wider(names_from = ticker,
#                      values_from = acao)
# 
# 
# ibov_filtrado = ibov %>% dplyr::mutate(month = month(date),
#                        year = year(date),
#                        day = day(date)) %>% 
#   dplyr::group_by(symbol) %>% 
#   dplyr::filter(all(year >= 2013)) %>%
#   dplyr::ungroup() %>% 
#   dplyr::select(-c(month, year, day)) %>% 
#   dplyr::mutate(date = make_date(year(date),
#                                  month(date),
#                                  day(date))) %>% 
#   dplyr::arrange(date) %>% 
#   dplyr::rename(acao = adjusted, ticker = symbol) %>% 
#   dplyr::select(ticker, date, acao) %>% print()
# 
# ibov_pivot_wider = ibov_filtrado %>% 
#   tidyr::pivot_wider(names_from = ticker,
#                      values_from = acao)

# Salvando os dataframes trabalhados

# write_rds(ativos_filtrado, file.path(getwd(), 'data', 'ativos_filtrado.rds'))
# 
# write_rds(ativos_modelagem, file.path(getwd(), 'data', 'ativos_modelagem.rds'))
# 
# write_rds(ibov_filtrado, file.path(getwd(), 'data', 'ibov_filtrado.rds'))
# 
# write_rds(ibov_pivot_wider, file.path(getwd(), 'data', 'ibov_pivot_wider.rds'))

Para poupar tempo, após todo o processo de baixar os dados dos ativos e do Ibovespa através do “tq_get()”, pelo yahoo finance, salvamos os dataframes em arquivos rds.

ativos_filtrado = read_rds(file.path(path, 'data','ativos_filtrado.rds'))

ativos_modelagem = read_rds(file.path(path, 'data', 'ativos_modelagem.rds'))

ibov_filtrado = read_rds(file.path(path, 'data', 'ibov_filtrado.rds'))

ibov_pivot_wider = read_rds(file.path(path, 'data','ibov_pivot_wider.rds'))


Estes são os ativos que estarão sendo operados em nosso portfolio.

unique(ativos_filtrado$ticker)
##   [1] "PETR3.SA"  "ITUB3.SA"  "VALE3.SA"  "ABEV3.SA"  "WEGE3.SA"  "BBAS3.SA" 
##   [7] "BBDC3.SA"  "ITSA3.SA"  "ELET3.SA"  "VIVT3.SA"  "SUZB3.SA"  "B3SA3.SA" 
##  [13] "SBSP3.SA"  "RENT3.SA"  "JBSS3.SA"  "RADL3.SA"  "PRIO3.SA"  "TIMS3.SA" 
##  [19] "CPFE3.SA"  "GGBR3.SA"  "EQTL3.SA"  "EGIE3.SA"  "CMIG3.SA"  "UGPA3.SA" 
##  [25] "BRFS3.SA"  "CSAN3.SA"  "CPLE3.SA"  "CCRO3.SA"  "EMBR3.SA"  "ENEV3.SA" 
##  [31] "PSSA3.SA"  "TRPL3.SA"  "CSNA3.SA"  "HYPE3.SA"  "EQPA3.SA"  "BRKM3.SA" 
##  [37] "ENMT3.SA"  "TOTS3.SA"  "CGAS3.SA"  "LREN3.SA"  "CIEL3.SA"  "REDE3.SA" 
##  [43] "MULT3.SA"  "USIM3.SA"  "MDIA3.SA"  "STBP3.SA"  "GOAU3.SA"  "MGLU3.SA" 
##  [49] "SMTO3.SA"  "MRFG3.SA"  "BNBR3.SA"  "EKTR3.SA"  "BRAP3.SA"  "SLCE3.SA" 
##  [55] "CSMG3.SA"  "CYRE3.SA"  "FLRY3.SA"  "POMO3.SA"  "WHRL3.SA"  "ENAT3.SA" 
##  [61] "CEEB3.SA"  "UNIP3.SA"  "ODPV3.SA"  "ALPA3.SA"  "DXCO3.SA"  "BAZA3.SA" 
##  [67] "ARZZ3.SA"  "GRND3.SA"  "ECOR3.SA"  "BRSR3.SA"  "EQMA3B.SA" "FRAS3.SA" 
##  [73] "MRSA3B.SA" "BBSE3.SA"  "RAIL3.SA"  "CRFB3.SA"  "VBBR3.SA"  "HAPV3.SA" 
##  [79] "NEOE3.SA"  "VIVA3.SA"  "NTCO3.SA"  "SOMA3.SA"  "CURY3.SA"  "GMAT3.SA" 
##  [85] "RRRP3.SA"  "RDOR3.SA"  "VAMO3.SA"  "INTB3.SA"  "CMIN3.SA"  "ASAI3.SA" 
##  [91] "AESB3.SA"  "GGPS3.SA"  "CXSE3.SA"  "RECV3.SA"  "SMFT3.SA"  "TTEN3.SA" 
##  [97] "PORT3.SA"  "SRNA3.SA"  "AURE3.SA"  "ALOS3.SA"


Preparação dos Dados

ativos_lista = ativos_filtrado %>% 
  split(f = .$ticker) %>% 
  keep(function(x){
    x$date %>% head(1) == "2013-01-02"
  }) 

ativos_lista = ativos_lista %>% 
  map(function(x){
    x %>% 
      mutate(SMA7 = TTR::SMA(acao, n = 7),
             SMA15 = TTR::SMA(acao, n = 15)) %>% 
      na.omit() %>% 
      mutate(posicao = case_when(lag(SMA7) > lag(SMA15) ~ 1,
                                 lag(SMA7) < lag(SMA15) ~ -1,
                                 lag(SMA7) == lag(SMA15) ~ 0), 
             retorno = case_when(posicao == 1 ~ acao/lag(acao) - 1,
                                 posicao == 0 ~ 0,
                                 posicao == -1 ~ lag(acao)/acao - 1)) %>% 
      na.omit()
  })


#Definindo funções:


Definindo a funcmax()

A função funcmax recebe dois parâmetros n1 e n2, que são usados para calcular médias móveis simples (SMA) de uma série temporal de preços de ações. A função realiza as seguintes operações:

Criação de uma lista de data frames modificados:

  • ativos_lista é uma lista de data frames onde cada data frame contém dados de uma ação.
  • Para cada data frame na lista ativos_lista, a função aplica uma série de transformações:
    • Calcula a média móvel simples de 7 e 15 períodos (ou qualquer valor de n1 e n2 fornecidos) para a coluna acao usando a função TTR::SMA.
    • Remove linhas com valores NA gerados pelo cálculo das SMAs com na.omit.
    • Cria uma nova coluna posicao que indica a posição da média móvel de curto prazo (SMA7) em relação à média móvel de longo prazo (SMA15):
      • 1 se SMA7 estiver acima de SMA15 (indica uma posição comprada).
      • -1 se SMA7 estiver abaixo de SMA15 (indica uma posição vendida).
      • 0 se SMA7 for igual a SMA15 (indica nenhuma posição).
    • Calcula a coluna retorno com base na posição:
      • Se posicao é 1, o retorno é calculado como a variação percentual positiva da ação.
      • Se posicao é -1, o retorno é calculado como a variação percentual negativa da ação.
      • Se posicao é 0, o retorno é 0.
    • Remove novamente linhas com valores NA gerados pelo cálculo das novas colunas com na.omit.

Criação do data frame da carteira:

  • Combina todos os data frames da lista em um único data frame.
  • Seleciona as colunas ticker, date e retorno.
  • Converte o formato do data frame de longo para largo, onde cada coluna representa os retornos de um ticker específico.
  • Adiciona uma nova coluna portfolio que calcula a média dos retornos de todas as ações para cada data (excluindo a coluna de datas) usando rowMeans.

Cálculo do retorno acumulado do portfólio:

  • Calcula o retorno acumulado do portfólio multiplicando os retornos diários incrementados em 1 e subtraindo 1 no final para obter o retorno total do portfólio. 
funcmax <- function(n1, n2){
  lista = ativos_lista %>% 
    map(function(x){
      x %>% 
        mutate(SMA7 = TTR::SMA(acao, n = n1),
               SMA15 = TTR::SMA(acao, n = n2)) %>% 
        na.omit() %>% 
        mutate(posicao = case_when(lag(SMA7) > lag(SMA15) ~ 1,
                                   lag(SMA7) < lag(SMA15) ~ -1,
                                   lag(SMA7) == lag(SMA15) ~ 0), 
               retorno = case_when(posicao == 1 ~ acao/lag(acao) - 1,
                                   posicao == 0 ~ 0,
                                   posicao == -1 ~ lag(acao)/acao - 1)) %>% 
        na.omit()
    })
  
  carteira <- lista %>% 
    bind_rows() %>% 
    select(ticker, date, retorno) %>% 
    pivot_wider(names_from = "ticker", values_from = "retorno") %>% 
    mutate(portfolio = rowMeans(.[-1], na.rm = TRUE))
  
  prod(carteira$portfolio + 1) - 1
}


Aqui definimos um dataframe contendo diversas combinações de médias móveis curtas e longas.

df_janelas_pre <- data.frame(sma_menor = sample(10:50, size = 10^6, replace = TRUE),
                         sma_maior = sample(10:50, size = 10^6, replace = TRUE)) %>%
filter(sma_menor < sma_maior)

head(df_janelas_pre)
##   sma_menor sma_maior
## 1        15        33
## 2        14        16
## 3        21        46
## 4        25        41
## 5        25        43
## 6        24        28


Backtest


Separando bases de treino e teste

Como temos ados de 2013 até o final de 2023, definimos que seria usado para treino os dados referentes a antes de 2020 e para teste os dados de depois de 2020.

ativos_treino = ativos_lista %>% 
  map(~.x %>%  filter(date < "2020-01-01"))

ativos_teste = ativos_lista %>% 
  map(~.x %>%  filter(date >= "2020-01-01"))


Treinamento do modelo

# df_janelas <- df_janelas %>% 
#   mutate(retorno = mapply(funcmaxTREINO, sma_menor, sma_maior))
# 
# df_janelas <- df_janelas %>% 
#   mutate(retornoanual=((retorno+1)^(1/6)-1))
# 
# df_janelas$retorno =  mapply(funcmax, df_janelas$sma_menor, df_janelas$sma_maior)
# 
# write_rds(df_janelas, file.path(path, "data", "df_janelas.rds"))

Novamente, para poupar tempo criamos arquivos rds com esses dfs gerados

df_janelas <- read_rds(file.path(path, "data", "df_janelas.rds"))
head(df_janelas)
##   sma_menor sma_maior  retorno retornoanual
## 1        43        47 1.580520    0.1084051
## 2        21        32 5.347833    0.2153398
## 3        28        50 3.649428    0.1845266
## 4        44        47 1.722562    0.1183525
## 5        36        50 2.507521    0.1430826
## 6        27        48 3.938750    0.1798360
df_janelas_max <- df_janelas %>% 
  filter(retorno == max(retorno))

df_janelas_max
##   sma_menor sma_maior  retorno retornoanual
## 1        22        30 5.832121    0.2215178

Podemos observar que o maior retorno foi obtido com as médias móveis nos períodos de 22 e 30 dias.


Teste do modelo

retteste<-funcmax(df_janelas_max$sma_menor[1],df_janelas_max$sma_maior[1])

retornoanualteste=((retteste+1)^(1/3)-1)

retornoanualteste
## [1] 0.8975148
funcretorno<-function(n1,n2){
  
  lista = ativos_teste %>% 
    map(function(x){
      x %>% 
        mutate(SMA7 = TTR::SMA(acao, n = n1),
               SMA15 = TTR::SMA(acao, n = n2)) %>% 
        na.omit() %>% 
        mutate(posicao = case_when(lag(SMA7) > lag(SMA15) ~ 1,
                                   lag(SMA7) < lag(SMA15) ~ -1,
                                   lag(SMA7) == lag(SMA15) ~ 0), 
               retorno = case_when(posicao == 1 ~ acao/lag(acao) - 1,
                                   posicao == 0 ~ 0,
                                   posicao == -1 ~ lag(acao)/acao - 1)) %>% 
        na.omit()
    })
  
  carteira<-lista %>% 
    bind_rows() %>% 
    select(ticker, date, retorno) %>% 
    pivot_wider(names_from = "ticker", values_from = "retorno") %>% 
    mutate(portfolio=rowMeans(.[-1],na.rm=T))
  
  return(carteira %>% select(date,portfolio))
  
}

df_backtest<-funcretorno(df_janelas_max$sma_menor[1],df_janelas_max$sma_maior[1])

head(df_backtest)
## # A tibble: 6 × 2
##   date       portfolio
##   <date>         <dbl>
## 1 2020-02-13  0.00125 
## 2 2020-02-14 -0.000497
## 3 2020-02-17 -0.000740
## 4 2020-02-18 -0.000839
## 5 2020-02-19  0.00172 
## 6 2020-02-20  0.00427


Análise de Retornos

df_backtest %>% 
  as.xts() %>% 
  chart.CumReturns()


Return.annualized(df_backtest %>% as.xts())
##                   portfolio
## Annualized Return 0.1943922


Importando dados do CDI e Ibov para fazer estatísticas e gráficos comparativos

# CDI
cdi<-read.csv2(file.path(path, 'data', 'cdi.csv'), stringsAsFactors = F)

cdi$Data <- dmy(cdi$Data)

cdi_filtrado = cdi %>% 
  filter(Data >= "2020-02-13") %>% 
  `colnames<-`(c('date', 'cdi')) %>% 
  mutate(cdi = as.numeric(gsub(",", ".", cdi))) %>% 
  filter(date <= "2023-12-28") %>% 
  mutate(cdi = cdi*0.01)

head(cdi_filtrado)
##         date        cdi
## 1 2020-02-13 0.00016137
## 2 2020-02-14 0.00016137
## 3 2020-02-17 0.00016137
## 4 2020-02-18 0.00016137
## 5 2020-02-19 0.00016137
## 6 2020-02-20 0.00016137
# Ibov
ibov_retornos <- ibov_filtrado %>% 
  filter((date <= "2023-12-28") & (date >= "2020-02-12")) %>% 
  select(date, ibov = acao) %>% 
  arrange(date) %>% 
  mutate(ibov = (ibov / lag(ibov)) - 1) %>% 
  na.omit()

head(ibov_retornos)
## # A tibble: 6 × 2
##   date           ibov
##   <date>        <dbl>
## 1 2020-02-13 -0.00867
## 2 2020-02-14 -0.0111 
## 3 2020-02-17  0.00811
## 4 2020-02-18 -0.00288
## 5 2020-02-19  0.0134 
## 6 2020-02-20 -0.0166

Preparando dados para plotagem

df_backtest_mod <- df_backtest %>% 
  mutate(ticker = "portfolio") %>% 
  select(ticker,date, retorno = portfolio)

cdi_filtrado <- cdi_filtrado %>% 
  mutate(ticker = "cdi") %>% 
  select(ticker, date, retorno = cdi)

ibov_retornos <- ibov_retornos %>% 
  mutate(ticker = "ibov") %>% 
  select(ticker, date, retorno = ibov)

df_plot = bind_rows(df_backtest_mod, cdi_filtrado, ibov_retornos)

df_plot_acumulado <- df_plot %>% 
  group_by(ticker) %>% 
  mutate(retorno_acumulado = cumprod(1 + retorno) - 1) %>% 
  mutate(retorno = retorno*100) %>% 
  mutate(retorno_acumulado = retorno_acumulado*100)


Gráficos

ggplot(df_plot_acumulado, aes(x = date, y = retorno_acumulado, color = ticker)) +
  geom_hline(yintercept = 0, color = "black") +
  geom_line() +
  labs(title = "Retorno Acumulado Estratégia x Benchmarks",
       subtitle = "Estratégia baseada no uso de médias móveis",
       x = "Data",
       y = "Retorno Acumulado",
       color = "Ticker",
       caption = "Fonte: Yahoo Finance, B3, BCB") +
  scale_x_date(date_labels = "%y-%m", breaks = scales::pretty_breaks(9)) +
  scale_y_continuous(labels = scales::percent_format(scale = 1), limits = c(-100, 100)) +
  theme_light() +
  theme(legend.position = "top") +
  scale_color_manual(values = c("portfolio" = "#0066FF",
                                "ibov" = "#4ead31",
                                "cdi" = "#fc0367"))


Estatísticas

Retorno acumulado Estratégia x Benchmarks

df_plot_acumulado %>%
  group_by(ticker) %>%
  summarise(retorno_acumulado = paste0(round(last(retorno_acumulado),1), "%"))
## # A tibble: 3 × 2
##   ticker    retorno_acumulado
##   <chr>     <chr>            
## 1 cdi       35.6%            
## 2 ibov      15%              
## 3 portfolio 97.2%


Retorno anualizado

df_port_xts = df_plot_estat %>% 
  filter(ticker == "portfolio")

df_port_xts <- xts(df_port_xts[, c("retorno", "retorno_acumulado")], order.by = df_port_xts$date)

retorno_anualizado_portfolio <- Return.annualized(df_port_xts$retorno)
valor_ra_port <- retorno_anualizado_portfolio[1,1]

df_ibov_xts = df_plot_estat %>% 
  filter(ticker == "ibov")

df_ibov_xts <- xts(df_ibov_xts[, c("retorno", "retorno_acumulado")], order.by = df_ibov_xts$date)

retorno_anualizado_ibov <- Return.annualized(df_ibov_xts$retorno)
valor_ra_ibov <- retorno_anualizado_ibov[1,1]


df_cdi_xts = df_plot_estat %>% 
  filter(ticker == "cdi")

df_cdi_xts <- xts(df_cdi_xts[, c("retorno", "retorno_acumulado")], order.by = df_cdi_xts$date)

retorno_anualizado_cdi <- Return.annualized(df_cdi_xts$retorno)
valor_ra_cdi <- retorno_anualizado_cdi[1,1]

df_ra <- data.frame("retorno_anualizado" = c("portfolio", "ibov", "cdi"),
                    "valor" = c(retorno_anualizado_portfolio, retorno_anualizado_ibov, retorno_anualizado_cdi)) %>% 
  mutate(valor = paste0(round(valor*100, 1),"%"))

df_ra
##   retorno_anualizado valor
## 1          portfolio 19.4%
## 2               ibov  3.7%
## 3                cdi  8.2%


Retornos positivos

positive_returns <- df_port_xts$retorno > 0
percentage_positive_returns <- mean(positive_returns) * 100
percentage_positive_returns
## [1] 53.47871


Gráfico calendário

table_calendar_return <- table.CalendarReturns(df_port_xts$retorno, digits = 1, as.perc = TRUE, geometric = TRUE)
table_calendar_return
##       Jan Feb Mar  Apr  May Jun  Jul  Aug  Sep  Oct Nov  Dec retorno
## 2020   NA 0.1 2.2 -1.5  0.5 0.3 -0.2  0.8 -0.1  0.0 0.1  0.0     2.1
## 2021 -0.5 0.8 0.0 -0.2  0.3 0.3  0.3  0.5  0.1  0.5 0.5  0.3     2.9
## 2022 -0.7 0.2 0.4 -0.2 -0.1 0.5 -0.1 -0.8  0.4  0.1 0.1  0.0    -0.2
## 2023  0.5 0.0 0.5  0.2  0.1 0.5 -0.3  0.6 -0.2 -0.3 0.8 -0.1     2.4


Métricas de Risco

Ferramentas quantitativas utilizadas para medir e avaliar o nível de risco associado a um investimento.
Para calculá-las utilizamos o pacote PerformanceAnalytics


VaR (Value at Risc)

indica a perda maxima esperada com uma probabilidade de 95%, ou seja, tem 95% de chance de que a perda maxima em um dia nao exceda o valor x

PerformanceAnalytics::VaR(df_port_xts$retorno, p = 0.95, method = "gaussian")
##        retorno
## VaR -0.0167313


Volatility Skewness

indica a assimetria da distribuicao da volatilidade dos retornos.
Ou seja, se o valor for positivo indica que existe uma chance maior de, caso ocorra uma volatilidade, dela ser ser positiva do que negativa, excedendo a media histórica.

PerformanceAnalytics::VolatilitySkewness(df_port_xts$retorno)
##          [,1]
## [1,] 2.119391


Gráfico de Drawdown

Representa as maiores quedas em relação aos picos históricos de um investimento.

dados_grafico_dd <- chart.Drawdown(df_port_xts, plot.engine = "ggplot2")
dados_grafico_dd  <- dados_grafico_dd$data
dados_grafico_dd <- dados_grafico_dd %>% 
  filter(security == "retorno")
ggplot(data = dados_grafico_dd,
       mapping = aes(x  = date,
                     y = returns,
                     color = security)) +
  geom_line() +
  scale_x_date(date_labels = "%y-%m",
               breaks = scales::pretty_breaks(9)) +
  scale_y_continuous(labels = scales::percent_format()) +
  labs(title = "Gráfico de drawdown da Estratégia",
       subtitle = "Estratégia baseada no uso de médias móveis",
       x = NULL,
       y = "Maiores quedas",
       color = NULL,
       caption = "Fonte: Yahoo Finance, B3") +
  theme_light() +
  theme(legend.position = "none")


Margens para melhoria e limitações da estratégia

. Utilizando bases de dados mais atualizadas poderíamos operar em margens de tempo menores.
. Poderíamos melhorar o método para calcular a combinação de médias móveis com melhores retornos, atualizando ela conforme o tempo ou calculando as melhores combinações para cada ativo.
. Poderíamos operar com outros tipos de ativos de alta liquidez, por exemplo criptomoedas, contratos futuros ou dólar.


. Nossa estratégia não leva em conta custos de transação nas operações.
. Nossa estratégia não leva em consideração problemas de liquidez em certos períodos, ou seja, conseguimos comprar ou vender qualquer ativo a qualquer momento.
- Isso pode não funcionar para certos momentos da história, por exemplo, na pandemia. Nossa estratégia detectou tendencias de queda na bolsa nesses períodos e operou vendido em diversos ativos, obtendo retornos bastante positivos em comparação com o Ibovespa.


Com isso encerramos a nossa apresentação. Obrigado!