packages_ap <- c(
"tidyverse", "lubridate", "readxl", "zoo",
"broom", "knitr", "kableExtra", "gridExtra"
)
for (pkg in packages_ap) {
if (!require(pkg, character.only = TRUE, quietly = TRUE)) {
install.packages(pkg, quiet = TRUE)
library(pkg, character.only = TRUE)
}
}
# Pacotes carregadosLista 2 - Econometria II (EPGE)
Sumário Executivo
Esta lista aborda dois grandes tópicos em econometria financeira: modelos de precificação de ativos e modelagem de séries temporais com componentes sazonais.
Questão 1: Asset Pricing
Aplicamos cinco modelos distintos para estimar o equity premium dos 25 Portfolios Formed on Size and Book-to-Market (5×5) de Fama-French (1947Q1-2025Q1):
- CCAPM (Consumption-based): Utiliza crescimento do consumo como fator de risco (γ=2, β=0.99)
- CAPM: Modelo clássico de um fator usando excesso de retorno de mercado
- Fama-French 3 Fatores: Adiciona fatores SMB (tamanho) e HML (valor)
- PCA: Extração estatística de fatores via componentes principais
- Fama-MacBeth: Metodologia em duas etapas para testar precificação cross-sectional
Principais Resultados: - CCAPM apresenta equity premium muito baixo, ilustrando o “equity premium puzzle” - FF3 oferece melhor ajuste (R² médio superior) comparado ao CAPM - PCA identifica componentes principais que explicam a maior parte da variação dos retornos
Questão 2: Modelagem SARIMA do IPCA
Aplicamos modelos SARIMA ao índice IPCA (jan/2000 - jan/2025):
- Diagnóstico: Série não-estacionária em nível, estacionária após diferenciação
- Modelos estimados: SARIMA(1,1,0)(1,0,0)₁₂ e SARIMA(1,1,0)(0,0,1)₁₂
- Previsões: Rolling window out-of-sample (fev/2023 - jan/2025)
Principais Resultados: - Ambos modelos apresentam parâmetros significativos e resíduos bem comportados - Modelo com SMA(1) apresenta melhor desempenho preditivo (menor RMSE) - Sazonalidade claramente capturada pelos componentes sazonais
Questão 1 - Modelos de Asset Pricing
Enunciado
Os modelos de asset pricing são fundamentais na teoria financeira, pois descrevem como os preços dos ativos são determinados no mercado. Começando pela perspectiva de um investidor que decide quanto poupar e consumir, e qual portfólio de ativos manter, a equação de precificação mais básica emerge das condições de primeira ordem desse problema. As equações desses modelos sugerem que o preço de um ativo deve ser o valor presente dos fluxos de caixa futuros esperados, utilizando a utilidade marginal do investidor, levando em conta o risco e a preferência temporal.
Com base nesses conceitos, aplicaremos os modelos para estimar o equity premium (excesso de retorno) dos 25 Portfolios Formed on Size and Book-to-Market (5×5) de Fama e French (disponíveis em Kenneth R. French - Data Library) conforme solicitado em cada item a seguir, utilizando dados trimestrais para o período entre 1947Q1 a 2025Q1. Quanto aos dados referentes ao consumo, utilizaremos o consumo de não duráveis em bases reais per capita norte-americano. As séries necessárias podem ser baixadas no site do Bureau of Economic Analysis (BEA), na seção das NIPA Tables. As séries a serem utilizadas estão nas tabelas 2.1 (linha 40), 2.3.4 e 2.3.5.
Setup e Carregamento de Pacotes - Asset Pricing
Carregamento e Preparação dos Dados
Dados de Consumo do BEA
# Carregamento dos dados do BEA
# Tabela 2.3.5 - Consumo nominal
bea_235_raw <- read_excel("BEATable2_3_5.xlsx", sheet = 1, col_names = FALSE)
# Encontrar linha com "Nondurable goods" (linha 8)
data_row_nondur <- which(grepl("^8$|^8 ", as.character(bea_235_raw[[1]])))[1]
# Extrair cabeçalhos de ano e trimestre
# Linha 6 tem os anos, linha 7 tem os quarters
years_row <- as.character(unlist(bea_235_raw[6, ]))
quarters_row <- as.character(unlist(bea_235_raw[7, ]))
# Começar da coluna 3 (primeira coluna de dados)
years_data <- years_row[3:length(years_row)]
quarters_data <- quarters_row[3:length(quarters_row)]
# Os anos aparecem apenas na primeira coluna de cada ano (e depois NA)
# Precisamos preencher os NAs com o último ano válido
current_year <- NA
for (i in seq_along(years_data)) {
if (!is.na(years_data[i]) && years_data[i] != "") {
current_year <- years_data[i]
} else {
years_data[i] <- current_year
}
}
# Criar nomes no formato YYYYQQ
col_names <- paste0(years_data, quarters_data)
# Filtrar apenas datas válidas (formato 2023Q1, 2024Q2, etc)
valid_indices <- grepl("^[0-9]{4}Q[1-4]$", col_names)
col_names <- col_names[valid_indices]
# Extrair dados de Nondurable goods (começando da coluna 3)
nondur_row <- unlist(bea_235_raw[data_row_nondur, ])
nondur_values_all <- as.numeric(nondur_row[3:length(nondur_row)])
nondur_values <- nondur_values_all[valid_indices]
# Criar dataframe
consumption_nominal <- data.frame(
Date_str = col_names,
Consumption_Nominal = nondur_values,
stringsAsFactors = FALSE
) %>%
mutate(Date = as.yearqtr(gsub("Q", " Q", Date_str))) # Formato "2023 Q1"
# Tabela 2.3.4 - Price indexes
bea_234_raw <- read_excel("BEATable2_3_4.xlsx", sheet = 1, col_names = FALSE)
data_row_nondur_price <- which(grepl("^8$|^8 ", as.character(bea_234_raw[[1]])))[1]
# Extrair cabeçalhos (mesma lógica da tabela 2.3.5)
years_row_234 <- as.character(unlist(bea_234_raw[6, ]))
quarters_row_234 <- as.character(unlist(bea_234_raw[7, ]))
years_data_234 <- years_row_234[3:length(years_row_234)]
quarters_data_234 <- quarters_row_234[3:length(quarters_row_234)]
# Preencher anos
current_year_234 <- NA
for (i in seq_along(years_data_234)) {
if (!is.na(years_data_234[i]) && years_data_234[i] != "") {
current_year_234 <- years_data_234[i]
} else {
years_data_234[i] <- current_year_234
}
}
col_names_234 <- paste0(years_data_234, quarters_data_234)
valid_indices_234 <- grepl("^[0-9]{4}Q[1-4]$", col_names_234)
col_names_234 <- col_names_234[valid_indices_234]
# Extrair price index
price_row <- unlist(bea_234_raw[data_row_nondur_price, ])
price_values_all <- as.numeric(price_row[3:length(price_row)])
price_values <- price_values_all[valid_indices_234]
price_index <- data.frame(
Date_str = col_names_234,
Price_Index = price_values,
stringsAsFactors = FALSE
) %>%
mutate(Date = as.yearqtr(gsub("Q", " Q", Date_str))) # Formato "2023 Q1"
# Tabela 2.1 - População (linha 40 tem população em milhares)
bea_21_raw <- read_excel("BEATable2_1.xlsx", sheet = 1, col_names = FALSE)
# A linha 40 está rotulada como "40" na coluna 1, mas está na linha 50 do Excel
# (há linhas de cabeçalho antes)
data_row_pop <- which(grepl("^40$|^40 ", as.character(bea_21_raw[[1]])))[1]
# Usar mesma lógica de extração
pop_row <- unlist(bea_21_raw[data_row_pop, ])
pop_values_all <- as.numeric(pop_row[3:length(pop_row)])
pop_values <- pop_values_all[valid_indices] # Usar mesmos índices válidos
population <- data.frame(
Date_str = col_names,
Population = pop_values,
stringsAsFactors = FALSE
) %>%
mutate(Date = as.yearqtr(gsub("Q", " Q", Date_str)))
# Consumo real per capita CORRIGIDO
consumption_data <- consumption_nominal %>%
inner_join(price_index, by = "Date") %>%
inner_join(population, by = "Date") %>%
mutate(
Consumption_Real = (Consumption_Nominal / Price_Index) * 100,
# Consumo real per capita: bilhões / (milhares/1000) = bilhões / milhões de pessoas
# Resultado em MIL DÓLARES per capita
Consumption_Real_PC = Consumption_Real / (Population / 1000)
) %>%
filter(Date >= as.yearqtr("1947 Q1") & Date <= as.yearqtr("2025 Q1")) %>%
select(Quarter = Date, Consumption_Real_PC, Consumption_Real, Population, Consumption_Nominal, Price_Index)
# Dados carregados com sucessoVisualização dos Dados do BEA
Descrição das variáveis:
- Quarter: Trimestre no formato YYYY QX
- Consumption_Real_PC: Consumo real per capita (mil dólares)
- Consumption_Real: Consumo real agregado (bilhões de dólares)
- Population: População (milhares de pessoas)
- Consumption_Nominal: Consumo nominal (bilhões de dólares)
- Price_Index: Índice de preços (base 100)
| Quarter | Consumption_Real_PC | Consumption_Real | Population | Consumption_Nominal | Price_Index |
|---|---|---|---|---|---|
| 1947 Q1 | 3.51 | 501.74 | 143143 | 74.9 | 14.93 |
| 1947 Q2 | 3.55 | 510.15 | 143790 | 76.9 | 15.07 |
| 1947 Q3 | 3.55 | 512.19 | 144449 | 78.6 | 15.35 |
| 1947 Q4 | 3.49 | 506.01 | 145122 | 80.0 | 15.81 |
| 1948 Q1 | 3.48 | 507.16 | 145709 | 81.5 | 16.07 |
| 1948 Q2 | 3.50 | 512.63 | 146289 | 83.2 | 16.23 |
| 1948 Q3 | 3.47 | 509.52 | 146921 | 83.5 | 16.39 |
| 1948 Q4 | 3.50 | 516.09 | 147607 | 83.7 | 16.22 |
| 1949 Q1 | 3.49 | 517.81 | 148254 | 82.7 | 15.97 |
| 1949 Q2 | 3.49 | 519.47 | 148847 | 81.9 | 15.77 |
| Quarter | Consumption_Real_PC | Consumption_Real | Population | Consumption_Nominal | Price_Index | |
|---|---|---|---|---|---|---|
| 304 | 2022 Q4 | 9.85 | 3302.93 | 335382 | 3896.6 | 117.97 |
| 305 | 2023 Q1 | 9.87 | 3317.37 | 336018 | 3926.9 | 118.37 |
| 306 | 2023 Q2 | 9.87 | 3322.14 | 336685 | 3938.6 | 118.56 |
| 307 | 2023 Q3 | 9.92 | 3349.27 | 337500 | 4000.9 | 119.46 |
| 308 | 2023 Q4 | 9.98 | 3375.98 | 338360 | 4032.2 | 119.44 |
| 309 | 2024 Q1 | 9.92 | 3362.74 | 339119 | 4010.0 | 119.25 |
| 310 | 2024 Q2 | 10.00 | 3397.94 | 339929 | 4072.5 | 119.85 |
| 311 | 2024 Q3 | 10.09 | 3438.21 | 340637 | 4107.6 | 119.47 |
| 312 | 2024 Q4 | 10.17 | 3470.29 | 341164 | 4143.7 | 119.40 |
| 313 | 2025 Q1 | 10.21 | 3488.97 | 341590 | 4196.5 | 120.28 |
| Variavel | Media | Desvio_Padrao | Minimo | Maximo |
|---|---|---|---|---|
| Consumo Real Per Capita | 6.02 | 1.86 | 3.45 | 10.21 |
| População (milhares) | 244050.36 | 58221.98 | 143143.00 | 341590.00 |
| Índice de Preços | 56.30 | 33.34 | 14.93 | 120.28 |
Portfólios Fama-French 25 (5x5)
# Ler arquivo CSV dos 25 portfolios 5x5 Size and B/M
ff_lines <- readLines("25_Portfolios_5x5.csv")
# Encontrar início dos dados mensais
start_line <- which(grepl("^\\s*,SMALL LoBM", ff_lines))[1] + 1
# Encontrar onde os dados terminam (antes de "Annual" ou nova seção)
# O arquivo pode ter várias seções - queremos apenas a primeira
end_markers <- which(grepl(
"Annual|^\\s*$|Average Value Weighted Returns|Average Equal Weighted Returns",
ff_lines[start_line:length(ff_lines)]
))
if (length(end_markers) > 0) {
end_line <- start_line + end_markers[1] - 2
} else {
end_line <- length(ff_lines)
}
# Ler dados
ff_monthly <- read.csv(
text = paste(ff_lines[start_line:end_line], collapse = "\n"),
header = FALSE, stringsAsFactors = FALSE
)
# Remover linhas vazias e linhas com "Annual"
ff_monthly <- ff_monthly[!grepl("Annual|^\\s*$", ff_monthly[, 1]), ]
ff_monthly <- ff_monthly[nchar(trimws(ff_monthly[, 1])) == 6, ] # Apenas YYYYMM
# Definir nomes das colunas
colnames(ff_monthly) <- c("Date", paste0("P", 1:25))
# Converter data e valores
ff_monthly <- ff_monthly %>%
mutate(
Date = as.Date(paste0(Date, "01"), format = "%Y%m%d"),
across(starts_with("P"), as.numeric)
) %>%
filter(!is.na(Date))
# Converter para trimestral usando retornos compostos
ff_monthly <- ff_monthly %>%
mutate(
Year = year(Date),
Quarter = quarter(Date),
YearQuarter = as.yearqtr(Date)
)
# Agregar para trimestral: (1+r1)(1+r2)(1+r3) - 1
# Nota: os retornos já estão em % nos dados originais, então dividimos por 100,
# calculamos o produto, e multiplicamos por 100 no final UMA NICA VEZ
ff_quarterly <- ff_monthly %>%
group_by(YearQuarter) %>%
summarise(across(
starts_with("P"),
~ ifelse(all(!is.na(.)), (prod(1 + . / 100) - 1) * 100, NA_real_)
), .groups = "drop") %>%
rename(Quarter = YearQuarter) %>%
filter(Quarter >= as.yearqtr("1947 Q1") & Quarter <= as.yearqtr("2025 Q1")) %>%
ungroup()
cat("Portfólios trimestrais:", nrow(ff_quarterly), "trimestres - Período:", as.character(min(ff_quarterly$Quarter)), "a", as.character(max(ff_quarterly$Quarter)), "\n")Portfólios trimestrais: 313 trimestres - Período: 1947 Q1 a 2025 Q1 Visualização dos Portfólios Fama-French
Descrição dos dados:
- 25 Portfolios Formed on Size and Book-to-Market (5×5)
- Ordenação dupla: 5 quintis de Tamanho (Market Equity) × 5 quintis de Book-to-Market
- Total: 5 × 5 = 25 portfólios distintos
- P1 a P5: Pequenas empresas, do menor ao maior B/M
- P21 a P25: Grandes empresas, do menor ao maior B/M
- Retornos em % ao trimestre (agregação geométrica de retornos mensais)
| Quarter | P1 | P2 | P3 | P4 | P5 | P6 | P7 | P8 | P9 | P10 | P11 | P12 | P13 | P14 | P15 | P16 | P17 | P18 | P19 | P20 | P21 | P22 | P23 | P24 | P25 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1947 Q1 | 1.3638 | -0.9658 | -3.1092 | 1.1195 | 3.6076 | -0.5583 | -1.9981 | 1.0741 | -1.6543 | 2.2709 | -2.3088 | -1.5196 | -1.9293 | -1.7944 | 2.7868 | -4.2820 | -1.7482 | -3.1170 | 0.7109 | -0.6586 | 1.2858 | -2.6817 | -2.3166 | -1.5210 | -3.2998 |
| 1947 Q2 | -10.9157 | -13.3635 | -9.8950 | -11.4121 | -9.0724 | -11.5351 | -5.6077 | -7.2369 | -7.2973 | -6.6614 | -10.4268 | -6.3806 | -6.5168 | -3.3511 | -5.1395 | -3.1947 | -4.4120 | -3.8292 | -8.1806 | -4.9963 | 0.2733 | 1.6759 | -2.8638 | 5.0963 | -1.8145 |
| 1947 Q3 | 8.6637 | 4.5249 | 3.4809 | 3.5805 | 9.7037 | 1.4044 | 7.5874 | 4.7050 | 6.5653 | 9.6582 | 7.5714 | 1.0960 | 7.0423 | 6.6337 | 8.4467 | 1.5653 | 1.8732 | 7.1358 | 5.0979 | 8.8058 | 0.9785 | 1.0722 | -0.2537 | 2.2055 | 3.4635 |
| 1947 Q4 | -5.8113 | 2.2603 | 2.0012 | 2.1558 | 5.5023 | -3.3611 | 1.6979 | 1.8376 | -0.3112 | 6.5555 | 4.4536 | 3.2169 | 0.9851 | 3.8144 | 1.9757 | 5.0147 | 3.1222 | 3.3497 | 0.2853 | -0.6833 | 1.6860 | 6.2377 | 4.7775 | 4.8852 | 9.7100 |
| 1948 Q1 | 6.9392 | -0.1354 | -1.0110 | 1.1433 | 5.2215 | -5.3257 | -2.0915 | -1.9946 | 2.8333 | 4.6885 | -2.4905 | -4.5020 | 1.9262 | -0.2601 | 6.2373 | -0.0890 | -1.5664 | 1.6523 | 4.0640 | 5.9790 | -2.6236 | 0.7531 | -1.1777 | 0.8380 | 1.8653 |
Séries Temporais dos Portfólios Selecionados
| Portfolio | Media | Desvio_Padrao | Minimo | Maximo | |
|---|---|---|---|---|---|
| P1 | P1 | 2.1877 | 15.1961 | -39.8975 | 54.3350 |
| P2 | P2 | 3.3672 | 12.8892 | -33.0183 | 51.5552 |
| P3 | P3 | 3.3909 | 11.6129 | -28.6682 | 44.0648 |
| P4 | P4 | 3.9598 | 11.2028 | -36.6349 | 43.5258 |
| P5 | P5 | 4.5579 | 12.3972 | -40.2378 | 53.4766 |
Fatores Fama-French
# Verificar se arquivo existe na pasta FF_Factors_temp
if (file.exists("FF_Factors_temp/F-F_Research_Data_Factors.csv")) {
ff_factors_lines <- readLines("FF_Factors_temp/F-F_Research_Data_Factors.csv")
} else {
stop("Arquivo F-F_Research_Data_Factors.csv não encontrado na pasta FF_Factors_temp")
}
# Encontrar início dos dados mensais
start_line <- which(grepl("^\\s*[0-9]{6}", ff_factors_lines))[1]
# Encontrar linha "Annual" para parar
annual_line <- which(grepl("Annual", ff_factors_lines))[1]
if (!is.na(annual_line)) {
end_line <- annual_line - 1
} else {
end_line <- length(ff_factors_lines)
}
# Ler dados mensais
ff_factors_monthly <- read.csv(
text = paste(ff_factors_lines[start_line:end_line], collapse = "\n"),
stringsAsFactors = FALSE
)
# Processar - renomear colunas para evitar problemas
colnames(ff_factors_monthly) <- c("Date", "Mkt_RF", "SMB", "HML", "RF")
ff_factors_monthly <- ff_factors_monthly %>%
filter(nchar(trimws(Date)) == 6) %>%
mutate(
Date = as.Date(paste0(Date, "01"), format = "%Y%m%d"),
across(-Date, as.numeric),
YearQuarter = as.yearqtr(Date)
) %>%
filter(!is.na(Date))
# Converter para trimestral (retornos compostos)
ff_factors_quarterly <- ff_factors_monthly %>%
group_by(YearQuarter) %>%
summarise(
Mkt_RF = (prod(1 + Mkt_RF / 100) - 1) * 100,
SMB = (prod(1 + SMB / 100) - 1) * 100,
HML = (prod(1 + HML / 100) - 1) * 100,
RF = (prod(1 + RF / 100) - 1) * 100,
.groups = "drop"
) %>%
rename(`Mkt-RF` = Mkt_RF, Quarter = YearQuarter) %>%
filter(Quarter >= as.yearqtr("1947 Q1") & Quarter <= as.yearqtr("2025 Q1"))
cat("Fatores trimestrais:", nrow(ff_factors_quarterly), "trimestres - Período:", as.character(min(ff_factors_quarterly$Quarter)), "a", as.character(max(ff_factors_quarterly$Quarter)), "\n")Fatores trimestrais: 313 trimestres - Período: 1947 Q1 a 2025 Q1 Visualização dos Fatores Fama-French
Descrição dos fatores:
- Mkt-RF: Excesso de retorno do mercado (market return - risk-free rate)
- SMB (Small Minus Big): Prêmio pelo fator tamanho
- HML (High Minus Low): Prêmio pelo fator valor (book-to-market)
- RF: Taxa livre de risco (risk-free rate)
- Todos em % ao trimestre
| Quarter | Mkt-RF | SMB | HML | RF |
|---|---|---|---|---|
| 1947 Q1 | -1.5361 | 1.2688 | 0.0053 | 0.0900 |
| 1947 Q2 | -0.7168 | -7.4167 | 0.5252 | 0.0900 |
| 1947 Q3 | 1.7656 | 3.4541 | 4.3034 | 0.1200 |
| 1947 Q4 | 3.4538 | -3.6197 | 4.9654 | 0.2001 |
| 1948 Q1 | -0.6878 | 0.8981 | 6.0154 | 0.2302 |
| 1948 Q2 | 11.1148 | -2.5503 | 5.7171 | 0.2502 |
| 1948 Q3 | -7.6971 | -2.7566 | -1.3768 | 0.2101 |
| 1948 Q4 | -0.8050 | -4.8717 | -5.4740 | 0.1200 |
| 1949 Q1 | 1.2336 | 2.2833 | 1.6634 | 0.2903 |
| 1949 Q2 | -4.6991 | -2.5383 | -5.0758 | 0.2903 |
| Fator | Media | Desvio_Padrao | Minimo | Maximo |
|---|---|---|---|---|
| Mkt-RF | 2.0724 | 8.1293 | -26.5150 | 23.1818 |
| SMB | 0.2544 | 4.9899 | -12.8864 | 15.9333 |
| HML | 0.9486 | 5.6925 | -22.2767 | 27.4500 |
| RF | 0.9612 | 0.7529 | 0.0000 | 3.8073 |
Visualização Temporal dos Fatores
Matriz de Correlação dos Fatores
| Mkt-RF | SMB | HML | RF | |
|---|---|---|---|---|
| Mkt-RF | 1.000 | 0.408 | -0.190 | -0.138 |
| SMB | 0.408 | 1.000 | -0.054 | -0.017 |
| HML | -0.190 | -0.054 | 1.000 | 0.069 |
| RF | -0.138 | -0.017 | 0.069 | 1.000 |
Interpretação da Correlação:
- Mkt-RF vs SMB: Correlação mediana indica que o fator tamanho captura risco ortogonal ao mercado
- Mkt-RF vs HML: Correlação baixa/negativa mostra independência entre fatores de mercado e valor
- SMB vs HML: Correlação próxima de zero confirma que capturam diferentes dimensões de risco
Inflação (CPI)
# Converter price index para taxa de inflação trimestral
cpi_quarterly <- price_index %>%
arrange(Date) %>%
mutate(Inflation_Rate = (Price_Index / lag(Price_Index) - 1) * 100) %>%
filter(Date >= as.yearqtr("1947 Q1") & Date <= as.yearqtr("2025 Q1")) %>%
select(Quarter = Date, Inflation_Rate, Price_Index)Combinação dos Dados
# Converter Quarter para character para fazer join
ff_quarterly_join <- ff_quarterly %>%
mutate(Quarter_chr = as.character(Quarter)) %>%
select(-Quarter)
ff_factors_join <- ff_factors_quarterly %>%
mutate(Quarter_chr = as.character(Quarter)) %>%
select(-Quarter)
consumption_join <- consumption_data %>%
mutate(Quarter_chr = as.character(Quarter)) %>%
select(-Quarter)
# Joins
data_full <- ff_quarterly_join %>%
left_join(ff_factors_join, by = "Quarter_chr") %>%
left_join(consumption_join %>% select(Quarter_chr, Consumption_Real_PC), by = "Quarter_chr") %>%
filter(!is.na(`Mkt-RF`), !is.na(Consumption_Real_PC)) %>%
mutate(Quarter = as.yearqtr(Quarter_chr)) %>%
select(Quarter, everything(), -Quarter_chr) %>%
arrange(Quarter)
# Adicionar inflação para calcular retornos reais
data_full <- data_full %>%
left_join(cpi_quarterly %>% select(Quarter, Inflation_Rate), by = "Quarter") %>%
filter(!is.na(Inflation_Rate))
# Calcular retornos reais para cada portfólio
# Nota: O enunciado sugere usar CPI mensal com a fórmula:
# ((1+r1/100)*(1+r2/100)*(1+r3/100))/((1+pi1/100)*(1+pi2/100)*(1+pi3/100)) - 1
# Porém, estamos usando o Price Index trimestral do BEA (Tabela 2.3.6), que já é agregado,
# então usamos a fórmula equivalente: (1 + r_trimestral) / (1 + inflacao_trimestral) - 1
for (i in 1:25) {
port_name <- paste0("P", i)
real_name <- paste0("P", i, "_Real")
data_full[[real_name]] <- ((1 + data_full[[port_name]] / 100) / (1 + data_full$Inflation_Rate / 100) - 1) * 100
}
# Calcular crescimento do consumo
data_full <- data_full %>%
mutate(Consumption_Growth = (Consumption_Real_PC / lag(Consumption_Real_PC) - 1) * 100) %>%
filter(!is.na(Consumption_Growth))
# Dataset final combinado(b) CAPM - Modelo de Um Fator
Enunciado: No site do French, são disponibilizados os dados referentes aos 3 fatores utilizados pelos autores em seu artigo de 1993, a saber: o excesso de retorno de mercado, e os fatores HML e SMB (HML = High minus Low: diferença entre os retornos das firmas com alto book-to-market e baixo book-to-market; SMB = Small minus Big: diferença entre empresas pequenas e grandes). Obtenha a planilha com esses dados, e estime o modelo CAPM para cada um dos portfólios, utilizando a variável (Mkt−Rf), e compute a média dos respectivos prêmios de risco.
Nota: O enunciado original menciona “HMB” e “SML”, mas as siglas corretas usadas por Fama e French são HML (High Minus Low) e SMB (Small Minus Big).
Resposta:
Fundamentação Teórica
O Capital Asset Pricing Model (CAPM), desenvolvido por Sharpe (1964), Lintner (1965) e Mossin (1966), é um dos modelos mais utilizados em finanças. Ele estabelece uma relação linear entre o retorno esperado de um ativo e seu risco sistemático (beta de mercado).
Hipóteses do modelo:
- Investidores são avessos ao risco e maximizam utilidade esperada
- Mercados são eficientes e sem fricções
- Todos os investidores têm o mesmo horizonte de investimento
- Existe uma taxa livre de risco para emprestar e tomar emprestado
Equação do CAPM:
Utilizando os dados dos 3 fatores disponibilizados por Fama e French, estimamos o modelo CAPM para cada portfólio usando a variável (Mkt-RF):
\[R^e_i = \alpha_i + \beta_i (R_M - R_f) + \epsilon_i\]
capm_results <- data.frame(Portfolio = paste0("P", 1:25))
# Calcular prêmio de mercado médio
market_premium <- mean(data_full$`Mkt-RF`, na.rm = TRUE)
for (i in 1:25) {
port_name <- paste0("P", i)
reg_data <- data_full %>%
transmute(
Excess_Ret = !!sym(port_name) - RF,
Mkt_RF = `Mkt-RF`
) %>%
filter(complete.cases(Excess_Ret, Mkt_RF))
if (nrow(reg_data) > 2) {
capm_model <- lm(Excess_Ret ~ Mkt_RF, data = reg_data)
capm_results[i, "Alpha"] <- coef(capm_model)[1]
capm_results[i, "Beta_Mkt"] <- coef(capm_model)[2]
capm_results[i, "R2"] <- summary(capm_model)$r.squared
# EP previsto pelo CAPM: Beta_Mkt × Prêmio de Mercado
capm_results[i, "EP_CAPM"] <- coef(capm_model)[2] * market_premium
} else {
capm_results[i, c("Alpha", "Beta_Mkt", "R2", "EP_CAPM")] <- NA
}
}RESULTADOS CAPM:Prêmio de Mercado Médio: 2.093 %Equity Premium Médio (CAPM): 2.3589 %Beta de Mercado Médio: 1.1271 R² Médio: 0.7587 | Portfolio | Alpha | Beta_Mkt | R2 | EP_CAPM |
|---|---|---|---|---|
| P1 | -1.9426 | 1.5329 | 0.6684 | 3.2083 |
| P2 | -0.3239 | 1.3340 | 0.7060 | 2.7920 |
| P3 | -0.0138 | 1.1952 | 0.6981 | 2.5014 |
| P4 | 0.6914 | 1.1277 | 0.6680 | 2.3602 |
| P5 | 1.1077 | 1.2089 | 0.6237 | 2.5302 |
| P6 | -1.1996 | 1.4338 | 0.7763 | 3.0009 |
| P7 | 0.0100 | 1.2333 | 0.7887 | 2.5812 |
| P8 | 0.4243 | 1.1150 | 0.7786 | 2.3336 |
| P9 | 0.6643 | 1.0649 | 0.7222 | 2.2288 |
| P10 | 0.8676 | 1.1940 | 0.6862 | 2.4990 |
| P11 | -0.7532 | 1.3291 | 0.8223 | 2.7818 |
| P12 | 0.1903 | 1.1325 | 0.8455 | 2.3702 |
| P13 | 0.3650 | 1.0247 | 0.7982 | 2.1447 |
| P14 | 0.7034 | 1.0535 | 0.7638 | 2.2050 |
| P15 | 0.8695 | 1.0992 | 0.6399 | 2.3006 |
| P16 | -0.3319 | 1.2256 | 0.8666 | 2.5651 |
| P17 | 0.0096 | 1.0736 | 0.8703 | 2.2469 |
| P18 | 0.4798 | 1.0094 | 0.8283 | 2.1126 |
| P19 | 0.6174 | 1.0179 | 0.7801 | 2.1305 |
| P20 | 0.5620 | 1.1320 | 0.6835 | 2.3693 |
| P21 | 0.0166 | 1.0106 | 0.8922 | 2.1151 |
| P22 | 0.0660 | 0.8982 | 0.8911 | 1.8799 |
| P23 | 0.5031 | 0.8222 | 0.7952 | 1.7208 |
| P24 | 0.1361 | 0.9046 | 0.7281 | 1.8933 |
| P25 | 0.4157 | 1.0036 | 0.6451 | 2.1005 |
Interpretação dos Resultados CAPM
Medidas apresentadas:
- Alpha: Intercepto da regressão (retorno anormal não explicado pelo mercado)
- Beta_Mkt: Sensibilidade do retorno ao excesso de retorno do mercado
- R²: Poder explicativo do modelo
- EP_CAPM: Equity premium previsto = Beta_Mkt × Prêmio de Mercado
Interpretação do Beta de Mercado:
- Beta > 1: Ativo mais volátil que o mercado (amplifica movimentos)
- Beta = 1: Ativo move-se com o mercado
- Beta < 1: Ativo menos volátil que o mercado
- Beta médio dos portfólios: 1.1271 Qualidade do Ajuste:
O R² médio de 0.7587 indica que o fator de mercado explica 75.87 % da variação dos retornos, um ajuste consideravelmente melhor que o CCAPM, mas ainda deixando espaço para outros fatores explicativos (como SMB e HML no modelo FF3).(c) Modelo de Três Fatores Fama-French
Enunciado: Em seguida, estime o modelo de três fatores proposto por Fama e French para cada um dos portfólios. Nesse caso, será necessário incluir os fatores HML e SMB na equação de regressão (ver definição em (b)). Novamente, compute o equity premium médio para cada um dos portfólios.
Resposta:
Fundamentação Teórica
O Modelo de Três Fatores de Fama-French (1993) expande o CAPM adicionando dois fatores que capturam anomalias documentadas empiricamente:
1. Fator SMB (Small Minus Big):
- Capta o “efeito tamanho” - empresas pequenas historicamente apresentam retornos superiores
- SMB = Retorno médio de empresas pequenas - Retorno médio de empresas grandes
2. Fator HML (High Minus Low):
- Capta o “efeito valor” - empresas com alto book-to-market (value) tendem a superar growth
- HML = Retorno médio de empresas com alto B/M - Retorno médio de empresas com baixo B/M
Motivação: Fama e French argumentam que esses fatores representam riscos sistemáticos adicionais não capturados pelo CAPM, relacionados a distress financeiro (HML) e problemas de acesso a capital (SMB).
Estimamos o modelo de três fatores proposto por Fama e French para cada portfólio, incluindo os fatores HML e SMB:
\[R^e_i = \alpha_i + \beta_{M,i} (R_M - R_f) + \beta_{SMB,i} SMB + \beta_{HML,i} HML + \epsilon_i\]
ff3_results <- data.frame(Portfolio = paste0("P", 1:25))
# Calcular prêmios médios dos fatores
market_premium <- mean(data_full$`Mkt-RF`, na.rm = TRUE)
smb_premium <- mean(data_full$SMB, na.rm = TRUE)
hml_premium <- mean(data_full$HML, na.rm = TRUE)
for (i in 1:25) {
port_name <- paste0("P", i)
reg_data <- data_full %>%
transmute(
Excess_Ret = !!sym(port_name) - RF,
Mkt_RF = `Mkt-RF`,
SMB = SMB,
HML = HML
) %>%
filter(complete.cases(Excess_Ret, Mkt_RF, SMB, HML))
if (nrow(reg_data) > 4) {
ff3_model <- lm(Excess_Ret ~ Mkt_RF + SMB + HML, data = reg_data)
ff3_results[i, "Alpha"] <- coef(ff3_model)[1]
ff3_results[i, "Beta_Mkt"] <- coef(ff3_model)[2]
ff3_results[i, "Beta_SMB"] <- coef(ff3_model)[3]
ff3_results[i, "Beta_HML"] <- coef(ff3_model)[4]
ff3_results[i, "R2"] <- summary(ff3_model)$r.squared
# EP previsto pelo FF3: Beta_Mkt × λ_Mkt + Beta_SMB × λ_SMB + Beta_HML × λ_HML
ff3_results[i, "EP_FF3"] <- coef(ff3_model)[2] * market_premium +
coef(ff3_model)[3] * smb_premium +
coef(ff3_model)[4] * hml_premium
} else {
ff3_results[i, c("Alpha", "Beta_Mkt", "Beta_SMB", "Beta_HML", "R2", "EP_FF3")] <- NA
}
}RESULTADOS FAMA-FRENCH 3 FATORES:Prêmio de Mercado: 2.093 %Prêmio SMB: 0.2758 %Prêmio HML: 0.953 %Equity Premium Médio (FF3): 2.5565 %R² Médio: 0.9199 | Portfolio | Alpha | Beta_Mkt | Beta_SMB | Beta_HML | R2 | EP_FF3 |
|---|---|---|---|---|---|---|
| P1 | -1.2301 | 1.1119 | 1.5418 | -0.2693 | 0.8866 | 2.4958 |
| P2 | -0.0783 | 1.0100 | 1.3335 | 0.0678 | 0.9276 | 2.5464 |
| P3 | -0.0669 | 0.9298 | 1.2162 | 0.2866 | 0.9460 | 2.5546 |
| P4 | 0.3675 | 0.9042 | 1.1598 | 0.4949 | 0.9553 | 2.6841 |
| P5 | 0.4604 | 1.0094 | 1.2100 | 0.7673 | 0.9454 | 3.1776 |
| P6 | -0.5122 | 1.1125 | 1.1068 | -0.3360 | 0.9379 | 2.3134 |
| P7 | 0.1819 | 0.9919 | 0.9992 | 0.0607 | 0.9517 | 2.4094 |
| P8 | 0.2228 | 0.9555 | 0.8132 | 0.3263 | 0.9427 | 2.5351 |
| P9 | 0.1686 | 0.9402 | 0.8007 | 0.5623 | 0.9506 | 2.7245 |
| P10 | 0.1109 | 1.0763 | 0.8959 | 0.7931 | 0.9568 | 3.2557 |
| P11 | -0.0840 | 1.0829 | 0.7785 | -0.3869 | 0.9404 | 2.1126 |
| P12 | 0.2425 | 0.9835 | 0.6426 | 0.0863 | 0.9335 | 2.3180 |
| P13 | 0.0578 | 0.9328 | 0.5618 | 0.3617 | 0.9232 | 2.4519 |
| P14 | 0.1089 | 1.0025 | 0.5183 | 0.5859 | 0.9374 | 2.7995 |
| P15 | 0.0044 | 1.0318 | 0.7235 | 0.8465 | 0.9124 | 3.1657 |
| P16 | 0.2376 | 1.0738 | 0.4037 | -0.3811 | 0.9338 | 1.9956 |
| P17 | -0.0414 | 1.0010 | 0.3504 | 0.1115 | 0.9042 | 2.2979 |
| P18 | 0.0706 | 0.9982 | 0.2494 | 0.3819 | 0.9017 | 2.5219 |
| P19 | 0.0608 | 1.0019 | 0.3427 | 0.5200 | 0.9067 | 2.6871 |
| P20 | -0.3318 | 1.1283 | 0.4514 | 0.8153 | 0.8890 | 3.2630 |
| P21 | 0.3689 | 1.0244 | -0.2334 | -0.3324 | 0.9541 | 1.7628 |
| P22 | -0.0460 | 0.9416 | -0.1401 | 0.0628 | 0.8997 | 1.9919 |
| P23 | 0.1386 | 0.9035 | -0.1876 | 0.2581 | 0.8440 | 2.0853 |
| P24 | -0.6822 | 1.0331 | -0.1781 | 0.6280 | 0.9009 | 2.7116 |
| P25 | -0.5337 | 1.1276 | -0.0940 | 0.7512 | 0.8169 | 3.0499 |
Interpretação dos Resultados FF3
Medidas apresentadas:
- Alpha: Intercepto (retorno anormal não explicado pelos três fatores)
- Beta_Mkt: Exposição ao risco de mercado
- Beta_SMB: Exposição ao fator tamanho (positivo = bias para small caps)
- Beta_HML: Exposição ao fator valor (positivo = bias para value stocks)
- R²: Poder explicativo do modelo
- EP_FF3: Equity premium previsto = Σ(Beta_i × Lambda_i) para os 3 fatores
Prêmios de Risco dos Fatores:
- Prêmio de Mercado (Mkt-RF): 2.093 % ao trimestre- Prêmio SMB: 0.2758 % ao trimestre- Prêmio HML: 0.953 % ao trimestreMelhoria em Relação ao CAPM:
- R² médio FF3: 0.9199 - R² médio CAPM: 0.7587 - Ganho no R²: 16.13 pontos percentuaisO modelo FF3 explica 91.99 % da variação dos retornos, uma melhoria substancial em relação ao CAPM. Isto confirma que os fatores SMB e HML capturam dimensões importantes do risco que não são explicadas apenas pelo mercado.(d) Análise de Componentes Principais
Enunciado: Usando os dados de retornos reais dos 25 Portfolios Formed on Size and Book-to-Market (5×5) de Fama e French, calcule os respectivos componentes principais (fatores). Quantos fatores explicam 60% ou mais da variação da base? Guarde-os.
Nota: Para calcular os retornos reais trimestrais, utilize a fórmula:
\[\left(\frac{(1+r_1/100)(1+r_2/100)(1+r_3/100)}{(1+\pi_1/100)(1+\pi_2/100)(1+\pi_3/100)} - 1\right) \times 100\]
utilizando o Índice de Preços ao Consumidor (CPI) em % ao mês como indicador de inflação.
Resposta:
Fundamentação Teórica
Análise de Componentes Principais (PCA) é uma técnica estatística que transforma variáveis correlacionadas em um conjunto menor de variáveis não correlacionadas chamadas componentes principais.
Objetivo no contexto de Asset Pricing:
- Identificar fatores de risco comuns que explicam a variação dos retornos
- Reduzir a dimensionalidade (de 25 portfólios para poucos fatores)
- Os componentes principais são combinações lineares dos retornos originais
- O primeiro PC capta a maior variação, o segundo PC capta a segunda maior variação ortogonal ao primeiro, e assim por diante
Diferença entre PCA e FF3:
- FF3 usa fatores construídos economicamente (mercado, tamanho, valor)
- PCA usa fatores construídos estatisticamente (máxima variação)
- PCA é puramente data-driven, sem teoria econômica a priori
Usando os dados de retornos reais dos 25 portfólios de Fama e French, calculamos os componentes principais (fatores) para determinar quantos fatores explicam 60% ou mais da variação da base.
# Matriz de retornos REAIS
returns_real_matrix <- data_full %>%
select(ends_with("_Real")) %>%
na.omit()
# PCA requer dados centrados para capturar corretamente a variação
returns_centered <- scale(returns_real_matrix, center = TRUE, scale = FALSE)
cat("Dados preparados para PCA: Retornos reais dos 25 portfólios, centralizados (média = 0), Período:", nrow(returns_centered), "observações\n\n")Dados preparados para PCA: Retornos reais dos 25 portfólios, centralizados (média = 0), Período: 311 observações# PCA (center=FALSE porque já centralizamos manualmente, scale=TRUE para padronizar variâncias)
pca_result <- prcomp(returns_centered, center = FALSE, scale. = TRUE)
# Variância explicada
variance_explained <- summary(pca_result)$importance[2, ]
cumvar <- cumsum(variance_explained)
# Quantos PCs para 60%?
n_pcs_60 <- which(cumvar >= 0.60)[1]
cat("Número de PCs para ≥ 60% da variância:", n_pcs_60, "- Variância acumulada:", round(cumvar[n_pcs_60] * 100, 2), "%\n\n")Número de PCs para ≥ 60% da variância: 1 - Variância acumulada: 85.14 %var_table <- data.frame(
PC = 1:10,
Var_Individual = round(variance_explained[1:10] * 100, 2),
Var_Acumulada = round(cumvar[1:10] * 100, 2)
)
kable(var_table,
col.names = c("Componente", "Variância (%)", "Variância Acumulada (%)"),
caption = "Variância Explicada por Componentes Principais"
) %>%
kable_styling(bootstrap_options = c("striped", "hover"))| Componente | Variância (%) | Variância Acumulada (%) | |
|---|---|---|---|
| PC1 | 1 | 85.14 | 85.14 |
| PC2 | 2 | 4.43 | 89.57 |
| PC3 | 3 | 3.87 | 93.44 |
| PC4 | 4 | 1.19 | 94.63 |
| PC5 | 5 | 0.77 | 95.40 |
| PC6 | 6 | 0.55 | 95.95 |
| PC7 | 7 | 0.52 | 96.46 |
| PC8 | 8 | 0.39 | 96.85 |
| PC9 | 9 | 0.38 | 97.23 |
| PC10 | 10 | 0.33 | 97.56 |
# Extrair e guardar os fatores (scores dos PCs)
pc_scores <- as.data.frame(pca_result$x[, 1:n_pcs_60])
colnames(pc_scores) <- paste0("PC", 1:n_pcs_60)
# Estimar equity premium usando PCA
pca_results <- data.frame(Portfolio = paste0("P", 1:25))
for (i in 1:25) {
port_name <- paste0("P", i, "_Real")
pca_data <- data.frame(
Excess_Ret = returns_real_matrix[[port_name]],
pc_scores
)
if (nrow(pca_data) > n_pcs_60 + 5) {
formula_pca <- as.formula(paste("Excess_Ret ~", paste(colnames(pc_scores), collapse = " + ")))
pca_model <- lm(formula_pca, data = pca_data)
pca_results[i, "R2_PCA"] <- summary(pca_model)$r.squared
pca_results[i, "EP_PCA"] <- mean(pca_data$Excess_Ret, na.rm = TRUE)
} else {
pca_results[i, c("R2_PCA", "EP_PCA")] <- NA
}
}
cat("\nEquity Premium Médio (PCA):", round(mean(pca_results$EP_PCA, na.rm = TRUE), 4), "% - R² Médio (PCA):", round(mean(pca_results$R2_PCA, na.rm = TRUE), 4), "\n")
Equity Premium Médio (PCA): 2.8123 % - R² Médio (PCA): 0.8514 # Gráfico de variância
var_df <- data.frame(
PC = 1:min(15, length(variance_explained)),
Variance = variance_explained[1:min(15, length(variance_explained))],
Cumulative = cumvar[1:min(15, length(variance_explained))]
)
p1 <- ggplot(var_df, aes(x = PC)) +
geom_col(aes(y = Variance), fill = "#3498DB", alpha = 0.8) +
geom_line(aes(y = Cumulative), color = "#E74C3C", linewidth = 1.5) +
geom_point(aes(y = Cumulative), color = "#E74C3C", size = 3) +
geom_hline(yintercept = 0.60, linetype = "dashed", color = "#27AE60", linewidth = 1) +
geom_vline(xintercept = n_pcs_60, linetype = "dotted", color = "#E74C3C", linewidth = 1) +
annotate("text",
x = n_pcs_60 + 1, y = 0.5,
label = paste(n_pcs_60, "PCs\n≥60%"),
color = "#E74C3C", fontface = "bold", size = 4
) +
scale_y_continuous(labels = scales::percent, breaks = seq(0, 1, 0.1)) +
labs(
title = "Análise de Componentes Principais: Variância Explicada",
subtitle = "25 Portfolios Formed on Size and Book-to-Market (5×5)",
x = "Componente Principal",
y = "Proporção da Variância",
caption = "Barras azuis: variância individual | Linha vermelha: variância acumulada"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 11),
axis.title = element_text(face = "bold"),
panel.grid.minor = element_blank()
)
print(p1)
(e) Comparação dos Modelos
Enunciado: Compare os equity premiums estimados nos itens (a), (b), (c) e (d), bem como o ajuste dos modelos com base no R².
Resposta:
# Combinar resultados
comparison <- data.frame(
Portfolio = paste0("P", 1:25),
EP_CCAPM = ccapm_results$EP_SDF,
EP_CAPM = capm_results$EP_CAPM,
EP_FF3 = ff3_results$EP_FF3,
EP_PCA = pca_results$EP_PCA,
R2_CCAPM = ccapm_results$R2,
R2_CAPM = capm_results$R2,
R2_FF3 = ff3_results$R2,
R2_PCA = pca_results$R2_PCA
)
# Médias
summary_table <- data.frame(
Modelo = c("CCAPM", "CAPM", "FF3", "PCA"),
EP_Medio = c(
mean(ccapm_results$EP_SDF, na.rm = TRUE),
mean(capm_results$EP_CAPM, na.rm = TRUE),
mean(ff3_results$EP_FF3, na.rm = TRUE),
mean(pca_results$EP_PCA, na.rm = TRUE)
),
R2_Medio = c(
mean(ccapm_results$R2, na.rm = TRUE),
mean(capm_results$R2, na.rm = TRUE),
mean(ff3_results$R2, na.rm = TRUE),
mean(pca_results$R2_PCA, na.rm = TRUE)
)
)
kable(summary_table,
digits = 4,
col.names = c("Modelo", "Equity Premium Médio (%)", "R² Médio"),
caption = "Comparação entre Modelos"
) %>%
kable_styling(bootstrap_options = c("striped", "hover"))| Modelo | Equity Premium Médio (%) | R² Médio |
|---|---|---|
| CCAPM | 0.0231 | 0.0162 |
| CAPM | 2.3589 | 0.7587 |
| FF3 | 2.5565 | 0.9199 |
| PCA | 2.8123 | 0.8514 |
# Gráfico comparativo
comparison_long <- comparison %>%
select(Portfolio, EP_CCAPM, EP_CAPM, EP_FF3, EP_PCA) %>%
pivot_longer(-Portfolio, names_to = "Modelo", values_to = "EP") %>%
mutate(
Modelo = factor(Modelo,
levels = c("EP_CCAPM", "EP_CAPM", "EP_FF3", "EP_PCA"),
labels = c("CCAPM", "CAPM", "FF3", "PCA")
),
# Ordenar Portfolio como P1, P2, ..., P25
Portfolio = factor(Portfolio, levels = paste0("P", 1:25))
) %>%
filter(!is.na(EP))
p2 <- ggplot(comparison_long, aes(x = Portfolio, y = EP, color = Modelo, group = Modelo)) +
geom_line(linewidth = 1.2, alpha = 0.7) +
geom_point(size = 2.5, alpha = 0.9) +
scale_color_manual(values = c(
"CCAPM" = "#E74C3C",
"CAPM" = "#3498DB",
"FF3" = "#27AE60",
"PCA" = "#9B59B6"
)) +
labs(
title = "Equity Premium por Portfólio: Comparação de Modelos",
subtitle = "25 Portfólios Fama-French (1947-2025)",
x = "Portfólio",
y = "Equity Premium (%)",
color = "Modelo",
caption = "Ordenação: Small/Big × Low/High (ME × INV)"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 11),
axis.text.x = element_text(angle = 45, hjust = 1, size = 8),
axis.title = element_text(face = "bold"),
legend.position = "top",
panel.grid.minor = element_blank()
)
print(p2)
# Gráfico comparativo de R²
r2_comparison <- data.frame(
Portfolio = paste0("P", 1:25),
CCAPM = ccapm_results$R2,
CAPM = capm_results$R2,
FF3 = ff3_results$R2,
PCA = pca_results$R2_PCA
) %>%
pivot_longer(-Portfolio, names_to = "Modelo", values_to = "R2") %>%
mutate(
Modelo = factor(Modelo, levels = c("CCAPM", "CAPM", "FF3", "PCA")),
Portfolio = factor(Portfolio, levels = paste0("P", 1:25))
) %>%
filter(!is.na(R2))
p3 <- ggplot(r2_comparison, aes(x = Modelo, y = R2, fill = Modelo)) +
geom_boxplot(alpha = 0.7, outlier.shape = 21, outlier.size = 3) +
geom_jitter(width = 0.2, alpha = 0.3, size = 1.5) +
scale_fill_manual(values = c(
"CCAPM" = "#E74C3C",
"CAPM" = "#3498DB",
"FF3" = "#27AE60",
"PCA" = "#9B59B6"
)) +
scale_y_continuous(labels = scales::percent) +
labs(
title = "Comparação do Ajuste dos Modelos",
subtitle = "Distribuição do R² nos 25 Portfólios",
x = NULL,
y = "R² (Coeficiente de Determinação)",
caption = "Box plots mostram mediana, quartis e outliers | Pontos individuais mostram cada portfólio"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 11),
axis.title = element_text(face = "bold"),
legend.position = "none",
panel.grid.major.x = element_blank()
)
print(p3)
Interpretação:
CCAPM: Apresenta equity premium médio muito baixo, refletindo o conhecido “equity premium puzzle” - a teoria básica de consumo tende a subestimar substancialmente os retornos excedentes de mercado.
CAPM e FF3: Estimam equity premiums com base nos betas estimados e nos prêmios de risco dos fatores de mercado.
R² do CCAPM: Muito baixo, indicando que o crescimento do consumo explica pouca variação nos retornos dos portfólios.
R² do FF3: Significativamente maior, demonstrando que os fatores de Fama-French (mercado, tamanho e valor) capturam melhor a variação cross-sectional dos retornos.
(f) Fama-MacBeth
Enunciado: Agora, iremos replicar o trabalho de Fama e Macbeth na nossa amostra, o qual consistirá em duas etapas.
Etapa em painel: combine os dados das duas tabelas do Fama e French em um único dataframe e estruture os dados para a regressão. Em seguida, realize regressões dos excessos de retornos dos 25 portfólios em relação aos três fatores, estimando os coeficientes beta para cada fator e portfólio.
Etapa cross section: utilize os coeficientes beta obtidos como variáveis independentes para calcular estatísticas resumidas, como média e desvio padrão, e realizar testes de significância. Finalmente, faça a regressão dos retornos médios dos portfólios sobre os coeficientes beta estimados e discuta as implicações do modelo.
Resposta:
A metodologia de Fama-MacBeth (1973) consiste em um procedimento de duas etapas para testar se os fatores de risco são precificados no mercado.
Fundamentação Teórica
O procedimento Fama-MacBeth testa a relação linear entre risco (betas) e retorno esperado:
\[E[R_i - R_f] = \lambda_0 + \lambda_{Mkt} \beta_{i,Mkt} + \lambda_{SMB} \beta_{i,SMB} + \lambda_{HML} \beta_{i,HML}\]
Onde: - \(\lambda_j\) são os prêmios de risco dos fatores (o que queremos estimar) - \(\beta_{i,j}\) são as exposições do portfólio \(i\) ao fator \(j\)
Etapa 1: Regressões em Séries Temporais
Para cada portfólio \(i\), estimamos os betas através de regressões de séries temporais:
\[R_{i,t} - R_{f,t} = \alpha_i + \beta_{i,Mkt}(R_{M,t} - R_{f,t}) + \beta_{i,SMB} SMB_t + \beta_{i,HML} HML_t + \epsilon_{i,t}\]
# ETAPA 1: Estimar betas para cada portfólio via time-series regressions
# (Já estimados na seção FF3 - vamos reutilizar)
betas_fm <- ff3_results %>%
select(Portfolio, Beta_Mkt, Beta_SMB, Beta_HML)
# Resumo dos betas
beta_summary <- data.frame(
Fator = c("Beta_Mkt", "Beta_SMB", "Beta_HML"),
Média = c(
mean(betas_fm$Beta_Mkt, na.rm = TRUE),
mean(betas_fm$Beta_SMB, na.rm = TRUE),
mean(betas_fm$Beta_HML, na.rm = TRUE)
),
Desvio_Padrão = c(
sd(betas_fm$Beta_Mkt, na.rm = TRUE),
sd(betas_fm$Beta_SMB, na.rm = TRUE),
sd(betas_fm$Beta_HML, na.rm = TRUE)
),
Mínimo = c(
min(betas_fm$Beta_Mkt, na.rm = TRUE),
min(betas_fm$Beta_SMB, na.rm = TRUE),
min(betas_fm$Beta_HML, na.rm = TRUE)
),
Máximo = c(
max(betas_fm$Beta_Mkt, na.rm = TRUE),
max(betas_fm$Beta_SMB, na.rm = TRUE),
max(betas_fm$Beta_HML, na.rm = TRUE)
)
)
kable(beta_summary,
digits = 4,
caption = "Resumo dos Betas Estimados (Etapa 1)") %>%
kable_styling(bootstrap_options = c("striped", "hover"))| Fator | Média | Desvio_Padrão | Mínimo | Máximo |
|---|---|---|---|---|
| Beta_Mkt | 1.0123 | 0.0684 | 0.9035 | 1.1283 |
| Beta_SMB | 0.6106 | 0.5162 | -0.2334 | 1.5418 |
| Beta_HML | 0.2825 | 0.4027 | -0.3869 | 0.8465 |
Etapa 2: Regressões Cross-Sectional em Cada Período
Para cada período \(t\), regredimos os retornos excedentes dos 25 portfólios contra seus betas:
\[R_{i,t} - R_{f,t} = \lambda_{0,t} + \lambda_{Mkt,t} \beta_{i,Mkt} + \lambda_{SMB,t} \beta_{i,SMB} + \lambda_{HML,t} \beta_{i,HML} + \eta_{i,t}\]
Os prêmios de risco são as médias dos \(\lambda_t\) ao longo do tempo.
# ETAPA 2: Para cada período t, fazer regressão cross-sectional
# Preparar dados: retornos excedentes de cada portfólio em cada período
n_periods <- nrow(data_full)
n_portfolios <- 25
# Matriz de retornos excedentes (períodos × portfólios)
excess_returns_matrix <- matrix(NA, nrow = n_periods, ncol = n_portfolios)
for (i in 1:n_portfolios) {
port_name <- paste0("P", i)
if (port_name %in% names(data_full)) {
excess_returns_matrix[, i] <- data_full[[port_name]] - data_full$RF
}
}
# Matriz de betas (portfólios × fatores)
beta_matrix <- as.matrix(betas_fm[, c("Beta_Mkt", "Beta_SMB", "Beta_HML")])
# Para cada período, rodar regressão cross-sectional
lambdas_t <- matrix(NA, nrow = n_periods, ncol = 4)
colnames(lambdas_t) <- c("Lambda_0", "Lambda_Mkt", "Lambda_SMB", "Lambda_HML")
for (t in 1:n_periods) {
# Retornos excedentes no período t
y_t <- excess_returns_matrix[t, ]
# Remover NAs
valid_idx <- !is.na(y_t) & complete.cases(beta_matrix)
if (sum(valid_idx) > 4) { # Precisamos de observações suficientes
# Regressão cross-sectional: R_i,t - R_f,t ~ betas
cs_data <- data.frame(
Ret = y_t[valid_idx],
Beta_Mkt = beta_matrix[valid_idx, 1],
Beta_SMB = beta_matrix[valid_idx, 2],
Beta_HML = beta_matrix[valid_idx, 3]
)
cs_reg <- lm(Ret ~ Beta_Mkt + Beta_SMB + Beta_HML, data = cs_data)
lambdas_t[t, ] <- coef(cs_reg)
}
}
# Remover períodos com NA
lambdas_t_clean <- na.omit(lambdas_t)
# Prêmios de risco: média dos lambdas ao longo do tempo
lambda_means <- colMeans(lambdas_t_clean)
lambda_sd <- apply(lambdas_t_clean, 2, sd)
lambda_se <- lambda_sd / sqrt(nrow(lambdas_t_clean))
lambda_t_stat <- lambda_means / lambda_se
lambda_p_value <- 2 * (1 - pt(abs(lambda_t_stat), df = nrow(lambdas_t_clean) - 1))
# Tabela de resultados
fama_macbeth_results <- data.frame(
Fator = c("Intercepto (α)", "Market (Mkt-RF)", "Size (SMB)", "Value (HML)"),
Lambda = lambda_means,
Erro_Padrão = lambda_se,
T_Stat = lambda_t_stat,
P_Valor = lambda_p_value,
Significativo = ifelse(lambda_p_value < 0.01, "***",
ifelse(lambda_p_value < 0.05, "**",
ifelse(lambda_p_value < 0.10, "*", "")))
)
kable(fama_macbeth_results,
digits = 4,
caption = "Prêmios de Risco Estimados - Procedimento Fama-MacBeth",
col.names = c("Fator", "λ (% trimestral)", "Erro Padrão", "Estatística t", "P-valor", "Signif.")
) %>%
kable_styling(bootstrap_options = c("striped", "hover"))| Fator | λ (% trimestral) | Erro Padrão | Estatística t | P-valor | Signif. | |
|---|---|---|---|---|---|---|
| Lambda_0 | Intercepto (α) | 3.1482 | 0.8028 | 3.9216 | 0.0001 | *** |
| Lambda_Mkt | Market (Mkt-RF) | -1.0422 | 0.9167 | -1.1368 | 0.2565 | |
| Lambda_SMB | Size (SMB) | 0.2445 | 0.2919 | 0.8376 | 0.4029 | |
| Lambda_HML | Value (HML) | 0.9975 | 0.3330 | 2.9958 | 0.0030 | *** |
Interpretação dos Resultados
**Prêmios de Risco Estimados:**- Intercepto (α): 3.1482 % - Retorno anormal médio- Market Premium: -1.0422 % - Prêmio de risco de mercado- SMB Premium: 0.2445 % - Prêmio pelo fator tamanho- HML Premium: 0.9975 % - Prêmio pelo fator valor**Testes de Significância:**✗ Fator de mercado não é significativamente precificado (p = 0.2565 )✗ Fator SMB não é significativamente precificado (p = 0.4029 )✓ Fator HML é significativamente precificado (p < 0.003 )Visualização: Evolução dos Prêmios de Risco no Tempo
# Gráfico da evolução dos lambdas ao longo do tempo
lambdas_df <- data.frame(
Periodo = 1:nrow(lambdas_t_clean),
Lambda_Mkt = lambdas_t_clean[, 2],
Lambda_SMB = lambdas_t_clean[, 3],
Lambda_HML = lambdas_t_clean[, 4]
) %>%
pivot_longer(cols = starts_with("Lambda"), names_to = "Fator", values_to = "Lambda")
lambdas_df$Fator <- factor(lambdas_df$Fator,
levels = c("Lambda_Mkt", "Lambda_SMB", "Lambda_HML"),
labels = c("Market Premium", "SMB Premium", "HML Premium")
)
p_lambdas <- ggplot(lambdas_df, aes(x = Periodo, y = Lambda, color = Fator)) +
geom_line(alpha = 0.3) +
geom_smooth(se = TRUE, linewidth = 1.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "black") +
facet_wrap(~ Fator, ncol = 1, scales = "free_y") +
scale_color_manual(values = c("#E74C3C", "#3498DB", "#27AE60")) +
labs(
title = "Evolução dos Prêmios de Risco - Fama-MacBeth",
subtitle = "Lambdas estimados em cada período via regressões cross-sectional",
x = "Período (Trimestre)",
y = "Prêmio de Risco (% trimestral)",
caption = "Linhas suaves mostram tendência (LOESS). Área sombreada = intervalo de confiança 95%"
) +
theme_minimal() +
theme(
legend.position = "none",
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 11),
strip.text = element_text(face = "bold", size = 11)
)
print(p_lambdas)
Conclusão: O procedimento Fama-MacBeth permite testar se os fatores identificados nas regressões de séries temporais (Etapa 1) efetivamente explicam a variação cross-sectional dos retornos esperados (Etapa 2). Fatores com \(\lambda\) significativamente diferentes de zero são considerados precificados pelo mercado e representam fontes de risco sistemático que os investidores exigem compensação.
Questão 2 - Modelo SARIMA
Enunciado
Sabe-se que séries econômicas de preços ao consumidor, como o Índice Nacional de Preços ao Consumidor Amplo (IPCA), podem apresentar componentes sazonais em virtude de padrões de consumo recorrentes ao longo do ano. A presença desses efeitos sazonais torna apropriada a utilização de modelos autorregressivos integrados de médias móveis sazonais, conhecidos como SARIMA \((p,d,q)(P,D,Q)_{12}\), que permitem capturar simultaneamente a dinâmica de curto prazo e a estrutura sazonal da série.
Com base na série mensal do IPCA entre janeiro de 2000 e janeiro de 2025, responda aos itens a seguir.
Atenção: utilizem a série do número-índice, e não a da taxa de variação; ela pode ser facilmente encontrada no site do IPEADATA.
Setup e Carregamento de Pacotes
packages <- c(
"tidyverse", "lubridate", "knitr", "kableExtra",
"forecast", "tseries", "urca", "ggplot2", "gridExtra",
"jsonlite", "httr"
)
for (pkg in packages) {
if (!require(pkg, character.only = TRUE, quietly = TRUE)) {
install.packages(pkg, quiet = TRUE)
library(pkg, character.only = TRUE)
}
}
# Pacotes carregadosCarregamento dos Dados do IPCA
# Baixando dados do IPCA do IPEADATA
ipca_url <- "http://ipeadata.gov.br/api/odata4/ValoresSerie(SERCODIGO='PRECOS12_IPCA12')"
ipca_raw <- jsonlite::fromJSON(ipca_url)
ipca_data <- ipca_raw$value %>%
mutate(Date = as.Date(VALDATA)) %>%
select(Date, IPCA = VALVALOR) %>%
filter(Date >= as.Date("2000-01-01") & Date <= as.Date("2025-01-31")) %>%
arrange(Date)
ipca_ts <- ts(ipca_data$IPCA, start = c(2000, 1), frequency = 12)
# Série temporal criadaVisualização dos Dados do IPCA
Descrição da variável:
- IPCA: Índice Nacional de Preços ao Consumidor Amplo (número-índice)
- Fonte: IPEADATA
- Frequência: Mensal
- Nota: Utilizamos o número-índice, não a taxa de variação
| Date | IPCA |
|---|---|
| 2000-01-01 | 1598.41 |
| 2000-02-01 | 1600.49 |
| 2000-03-01 | 1604.01 |
| 2000-04-01 | 1610.75 |
| 2000-05-01 | 1610.91 |
| 2000-06-01 | 1614.62 |
| 2000-07-01 | 1640.62 |
| 2000-08-01 | 1662.11 |
| 2000-09-01 | 1665.93 |
| 2000-10-01 | 1668.26 |
| 2000-11-01 | 1673.60 |
| 2000-12-01 | 1683.47 |
| Date | IPCA | |
|---|---|---|
| 290 | 2024-02-01 | 6858.17 |
| 291 | 2024-03-01 | 6869.14 |
| 292 | 2024-04-01 | 6895.24 |
| 293 | 2024-05-01 | 6926.96 |
| 294 | 2024-06-01 | 6941.51 |
| 295 | 2024-07-01 | 6967.89 |
| 296 | 2024-08-01 | 6966.50 |
| 297 | 2024-09-01 | 6997.15 |
| 298 | 2024-10-01 | 7036.33 |
| 299 | 2024-11-01 | 7063.77 |
| 300 | 2024-12-01 | 7100.50 |
| 301 | 2025-01-01 | 7111.86 |
| Estatistica | Valor |
|---|---|
| Média | 3858.08 |
| Desvio Padrão | 1577.83 |
| Mínimo | 1598.41 |
| Máximo | 7111.86 |
| Mediana | 3497.70 |
(a) Estatísticas Descritivas e Testes de Raiz Unitária
Enunciado: Apresente a estatística descritiva da série, incluindo gráficos da série em nível e com as devidas transformações. Aplique os testes de raiz unitária de ADF e Phillips-Perron para verificar se a série é estacionária em nível, em log ou na primeira diferença do log. Além disso, plote a função de autocorrelação (FAC) e a função de autocorrelação parcial (FACP), para identificação das ordens p, q, P e Q do modelo.
Resposta:
Fundamentação Teórica sobre Estacionariedade
Série Estacionária: Uma série temporal é estacionária quando suas propriedades estatísticas (média, variância, autocorrelação) não mudam ao longo do tempo.
Por que importa?
- Modelos ARIMA/SARIMA requerem estacionariedade
- Séries não-estacionárias podem levar a “regressões espúrias”
- Preços geralmente são não-estacionários (têm tendência), mas retornos/variações são estacionários
Testes de Raiz Unitária:
- ADF (Augmented Dickey-Fuller):
- H₀: A série tem raiz unitária (não-estacionária)
- H₁: A série é estacionária
- P-valor < 0.05 → Rejeita H₀ → Série é estacionária
- Phillips-Perron:
- Similar ao ADF, mas robusto a autocorrelação e heterocedasticidade
- Mesma interpretação do p-valor
Transformações Comuns:
- Log: Estabiliza variância
- Primeira Diferença: Remove tendência linear
- Diferença Sazonal (lag=12): Remove sazonalidade
summary(ipca_ts) Min. 1st Qu. Median Mean 3rd Qu. Max.
1598 2574 3498 3858 5093 7112 par(mfrow = c(2, 2))
plot(ipca_ts, main = "IPCA - Nível", col = "#003366", lwd = 1.5)
grid()
log_ipca <- log(ipca_ts)
plot(log_ipca, main = "IPCA - Log", col = "#003366", lwd = 1.5)
grid()
diff_log_ipca <- diff(log_ipca)
plot(diff_log_ipca, main = "IPCA - Primeira Diferença", col = "#003366", lwd = 1.5)
grid()
abline(h = 0, col = "black", lty = 2)
diff12_log_ipca <- diff(log_ipca, lag = 12)
plot(diff12_log_ipca, main = "IPCA - Diferença Sazonal", col = "#003366", lwd = 1.5)
grid()
abline(h = 0, col = "black", lty = 2)
par(mfrow = c(1, 1))
# Testes de raiz unitária
adf_level <- adf.test(ipca_ts)
adf_log <- adf.test(log_ipca)
adf_diff <- adf.test(diff_log_ipca)
pp_level <- pp.test(ipca_ts)
pp_log <- pp.test(log_ipca)
pp_diff <- pp.test(diff_log_ipca)
par(mfrow = c(2, 2))
acf(diff_log_ipca, main = "FAC - Δlog(IPCA)", lag.max = 48, col = "#003366", lwd = 2)
pacf(diff_log_ipca, main = "FACP - Δlog(IPCA)", lag.max = 48, col = "#DC143C", lwd = 2)
diff_both <- diff(diff_log_ipca, lag = 12)
acf(diff_both, main = "FAC - Diferenças", lag.max = 48, na.action = na.pass, col = "#003366", lwd = 2)
pacf(diff_both, main = "FACP - Diferenças", lag.max = 48, na.action = na.pass, col = "#DC143C", lwd = 2)
par(mfrow = c(1, 1))Teste ADF:
- Nível - P-valor: 0.9846 - Log - P-valor: 0.4413 - Primeira Diferença do Log - P-valor: 0.01 Teste Phillips-Perron:
- Nível - P-valor: 0.99 - Log - P-valor: 0.6953 - Primeira Diferença do Log - P-valor: 0.01 Interpretação dos Testes e Gráficos
Testes de Raiz Unitária:
- Série em nível: NÃO é estacionária (p-valor = 0.9846 )
- Série em log: NÃO é estacionária (p-valor = 0.4413 )
- Primeira diferença do log: É estacionária (p-valor = 0.01 )Conclusão: A série IPCA em nível apresenta tendência estocástica (raiz unitária). Após aplicar a primeira diferença no logaritmo, a série torna-se estacionária, indicando que devemos usar d=1 no modelo SARIMA.
Funções de Autocorrelação (FAC e FACP):
- FAC: Mostra correlação entre a série e seus próprios lags
- FACP: Mostra correlação direta (removendo efeitos de lags intermediários)
- Picos em múltiplos de 12: Indicam forte componente sazonal
- Decaimento gradual na FAC: Sugere componente AR
- Corte abrupto na FACP: Ajuda a identificar ordem p
Identificação do Modelo:
Com base nos gráficos FAC/FACP, os modelos SARIMA(1,1,0)(1,0,0)₁₂ e SARIMA(1,1,0)(0,0,1)₁₂ são candidatos razoáveis, combinando:
- Componente AR de curto prazo (p=1)
- Diferenciação regular (d=1)
- Componente sazonal (P=1 ou Q=1 com período s=12)
(b) Estimação dos Modelos SARIMA
Enunciado: Estime, utilizando a amostra completa, os modelos SARIMA(1,1,0)(1,0,0)₁₂ e SARIMA(1,1,0)(0,0,1)₁₂. Para cada um deles, obtenha os parâmetros estimados, realize os diagnósticos (significância dos parâmetros, análise da ACF e PACF dos resíduos e verificação das raízes do polinômio característico) e compare os modelos em termos dos Critérios de Informação (AIC e BIC).
Resposta:
Modelo 1: SARIMA(1,1,0)(1,0,0)₁₂
# MODELO 1: SARIMA(1,1,0)(1,0,0)[12]
model1 <- Arima(log_ipca, order = c(1, 1, 0), seasonal = list(order = c(1, 0, 0), period = 12))
summary(model1)Series: log_ipca
ARIMA(1,1,0)(1,0,0)[12]
Coefficients:
ar1 sar1
0.8060 0.2270
s.e. 0.0376 0.0657
sigma^2 = 1.06e-05: log likelihood = 1295.33
AIC=-2584.65 AICc=-2584.57 BIC=-2573.54
Training set error measures:
ME RMSE MAE MPE MAPE
Training set 0.00077353 0.003238878 0.002341845 0.009619439 0.0286874
MASE ACF1
Training set 0.03888144 -0.1301346Parâmetros Estimados:
AR(1): 0.806 (SE: 0.0376 )SAR(1): 0.227 (SE: 0.0657 )Teste de Significância:
AR(1): t = 21.463 , p-valor = 0 SAR(1): t = 3.458 , p-valor = 5e-04 Raízes do Polinômio:
Raízes AR: 1.2407 Raízes SAR: 4.4052 Estacionário: TRUE # Extrair resíduos
res1 <- residuals(model1)
# 1. Gráfico de resíduos ao longo do tempo
p1 <- autoplot(res1, colour = "#003366") +
geom_hline(yintercept = 0, linetype = "dashed", color = "#DC143C", linewidth = 0.8) +
labs(
title = "Resíduos ao Longo do Tempo",
x = "Tempo",
y = "Resíduos"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 12),
panel.grid.minor = element_blank()
)
# 2. ACF dos resíduos
acf_data <- acf(res1, plot = FALSE, lag.max = 24)
acf_df <- data.frame(
Lag = acf_data$lag[-1],
ACF = acf_data$acf[-1]
)
ci <- qnorm((1 + 0.95) / 2) / sqrt(length(res1))
p2 <- ggplot(acf_df, aes(x = Lag, y = ACF)) +
geom_hline(yintercept = 0, color = "black", linewidth = 0.5) +
geom_hline(yintercept = c(-ci, ci), linetype = "dashed", color = "#DC143C", linewidth = 0.8) +
geom_segment(aes(xend = Lag, yend = 0), color = "#003366", linewidth = 1) +
geom_point(color = "#003366", size = 2) +
labs(
title = "ACF dos Resíduos",
x = "Lag",
y = "ACF"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 12),
panel.grid.minor = element_blank()
)
# 3. Histograma com curva normal
res_df <- data.frame(Residuos = as.numeric(res1))
p3 <- ggplot(res_df, aes(x = Residuos)) +
geom_histogram(aes(y = after_stat(density)), bins = 30, fill = "#003366", alpha = 0.7, color = "white") +
stat_function(
fun = dnorm,
args = list(mean = mean(res1, na.rm = TRUE), sd = sd(res1, na.rm = TRUE)),
color = "#DC143C",
linewidth = 1.2
) +
labs(
title = "Distribuição dos Resíduos",
x = "Resíduos",
y = "Densidade"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 12),
panel.grid.minor = element_blank()
)
# Combinar gráficos: 1 em cima, 2 embaixo
layout_matrix <- rbind(c(1, 1), c(2, 3))
grid.arrange(p1, p2, p3, layout_matrix = layout_matrix, heights = c(1, 1))
Teste de Ljung-Box:
P-valor = 0.0749 Resíduos não apresentam autocorrelação significativa (ruído branco).Modelo 2: SARIMA(1,1,0)(0,0,1)₁₂
# MODELO 2: SARIMA(1,1,0)(0,0,1)[12]
model2 <- Arima(log_ipca, order = c(1, 1, 0), seasonal = list(order = c(0, 0, 1), period = 12))
summary(model2)Series: log_ipca
ARIMA(1,1,0)(0,0,1)[12]
Coefficients:
ar1 sma1
0.8232 0.1775
s.e. 0.0336 0.0559
sigma^2 = 1.069e-05: log likelihood = 1294.07
AIC=-2582.14 AICc=-2582.05 BIC=-2571.02
Training set error measures:
ME RMSE MAE MPE MAPE
Training set 0.0007725418 0.003253154 0.002351729 0.009602919 0.02882324
MASE ACF1
Training set 0.03904555 -0.1386685Parâmetros Estimados:
AR(1): 0.8232 (SE: 0.0336 )SMA(1): 0.1775 (SE: 0.0559 )Teste de Significância:
AR(1): t = 24.498 , p-valor = 0 SMA(1): t = 3.172 , p-valor = 0.0015 Raízes do Polinômio:
Raízes AR: 1.2147 Raízes SMA: 5.6353 Estacionário e invertível: TRUE # Extrair resíduos
res2 <- residuals(model2)
# 1. Gráfico de resíduos ao longo do tempo
p1 <- autoplot(res2, colour = "#003366") +
geom_hline(yintercept = 0, linetype = "dashed", color = "#DC143C", linewidth = 0.8) +
labs(
title = "Resíduos ao Longo do Tempo",
x = "Tempo",
y = "Resíduos"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 12),
panel.grid.minor = element_blank()
)
# 2. ACF dos resíduos
acf_data <- acf(res2, plot = FALSE, lag.max = 24)
acf_df <- data.frame(
Lag = acf_data$lag[-1],
ACF = acf_data$acf[-1]
)
ci <- qnorm((1 + 0.95) / 2) / sqrt(length(res2))
p2 <- ggplot(acf_df, aes(x = Lag, y = ACF)) +
geom_hline(yintercept = 0, color = "black", linewidth = 0.5) +
geom_hline(yintercept = c(-ci, ci), linetype = "dashed", color = "#DC143C", linewidth = 0.8) +
geom_segment(aes(xend = Lag, yend = 0), color = "#003366", linewidth = 1) +
geom_point(color = "#003366", size = 2) +
labs(
title = "ACF dos Resíduos",
x = "Lag",
y = "ACF"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 12),
panel.grid.minor = element_blank()
)
# 3. Histograma com curva normal
res_df <- data.frame(Residuos = as.numeric(res2))
p3 <- ggplot(res_df, aes(x = Residuos)) +
geom_histogram(aes(y = after_stat(density)), bins = 30, fill = "#003366", alpha = 0.7, color = "white") +
stat_function(
fun = dnorm,
args = list(mean = mean(res2, na.rm = TRUE), sd = sd(res2, na.rm = TRUE)),
color = "#DC143C",
linewidth = 1.2
) +
labs(
title = "Distribuição dos Resíduos",
x = "Resíduos",
y = "Densidade"
) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 12),
panel.grid.minor = element_blank()
)
# Combinar gráficos: 1 em cima, 2 embaixo
layout_matrix <- rbind(c(1, 1), c(2, 3))
grid.arrange(p1, p2, p3, layout_matrix = layout_matrix, heights = c(1, 1))
Teste de Ljung-Box:
P-valor = 0.0634 Resíduos não apresentam autocorrelação significativa (ruído branco).Comparação dos Modelos
| Modelo | AIC | BIC | LogLik |
|---|---|---|---|
| SARIMA(1,1,0)(1,0,0)[12] | -2584.653 | -2573.541 | 1295.326 |
| SARIMA(1,1,0)(0,0,1)[12] | -2582.135 | -2571.024 | 1294.068 |
Interpretação:
- AIC prefere o Modelo 1 (SAR)- BIC prefere o Modelo 1 (SAR)
Critérios menores indicam melhor ajuste. BIC penaliza mais a complexidade do modelo.(c) Previsões Fora da Amostra
Enunciado: Considere agora apenas os dados até janeiro de 2023. Utilizando as duas especificações da questão anterior, obtenha previsões para o período de fevereiro de 2023 até janeiro de 2025. Neste exercício, utilize a previsão passo a passo: estime o modelo e obtenha a previsão para h=1; em seguida, incorpore a observação seguinte à amostra, reestime o modelo e gere a previsão para h=1, repetindo esse procedimento até o final do horizonte. Plote as previsões resultantes junto com os dados observados e compare o desempenho de cada modelo com base no RMSE.
Resposta:
# Previsões Rolling Forecast
ipca_train <- window(log_ipca, end = c(2023, 1))
ipca_test <- window(log_ipca, start = c(2023, 2))
rolling_forecast <- function(train_data, test_data, model_order, seasonal_order) {
forecasts <- numeric(length(test_data))
for (i in 1:length(test_data)) {
if (i == 1) {
current_data <- train_data
} else {
current_data <- ts(c(train_data, test_data[1:(i - 1)]),
start = start(train_data), frequency = 12
)
}
model <- Arima(current_data,
order = model_order,
seasonal = list(order = seasonal_order, period = 12)
)
forecast_h1 <- forecast(model, h = 1)
forecasts[i] <- as.numeric(forecast_h1$mean)
}
return(forecasts)
}
# Gerar previsões
forecasts_m1 <- rolling_forecast(ipca_train, ipca_test, c(1, 1, 0), c(1, 0, 0))
forecasts_m2 <- rolling_forecast(ipca_train, ipca_test, c(1, 1, 0), c(0, 0, 1))
rmse_m1 <- sqrt(mean((ipca_test - forecasts_m1)^2))
rmse_m2 <- sqrt(mean((ipca_test - forecasts_m2)^2))
# Métricas de desempenho
performance <- data.frame(
Modelo = c("SARIMA(1,1,0)(1,0,0)[12]", "SARIMA(1,1,0)(0,0,1)[12]"),
RMSE = c(rmse_m1, rmse_m2)
)
kable(performance, digits = 6) %>% kable_styling()| Modelo | RMSE |
|---|---|
| SARIMA(1,1,0)(1,0,0)[12] | 0.002595 |
| SARIMA(1,1,0)(0,0,1)[12] | 0.002565 |
if (rmse_m1 < rmse_m2) {
cat("\n✓ Modelo 1 apresenta melhor desempenho preditivo\n")
} else {
cat("\n✓ Modelo 2 apresenta melhor desempenho preditivo\n")
}
✓ Modelo 2 apresenta melhor desempenho preditivo# Análise de similaridade entre as previsões
correlation <- cor(forecasts_m1, forecasts_m2)
diff_forecasts <- forecasts_m1 - forecasts_m2
similarity_metrics <- data.frame(
Métrica = c(
"Correlação entre modelos",
"Diferença média (log)",
"Diferença máxima absoluta (log)",
"Desvio padrão das diferenças"
),
Valor = c(
round(correlation, 4),
round(mean(diff_forecasts), 6),
round(max(abs(diff_forecasts)), 6),
round(sd(diff_forecasts), 6)
)
)
cat("\n")kable(similarity_metrics,
caption = "Análise de Similaridade das Previsões",
col.names = c("Métrica", "Valor")
) %>%
kable_styling(bootstrap_options = c("striped", "hover"))| Métrica | Valor |
|---|---|
| Correlação entre modelos | 0.999900 |
| Diferença média (log) | 0.000017 |
| Diferença máxima absoluta (log) | 0.000511 |
| Desvio padrão das diferenças | 0.000254 |
Interpretação: As previsões são muito similares porque ambos os modelos compartilham a mesma estrutura não-sazonal (AR(1) + I(1)) e diferem apenas na componente sazonal (SAR(1) vs SMA(1)). Em previsões de um passo à frente (h=1), essas especificações sazonais tendem a produzir resultados próximos, especialmente quando o componente AR(1) domina a dinâmica de curto prazo.
Conclusões
Questão 1 - Asset Pricing
Principais Resultados:
- CCAPM (Consumption-based):
- O modelo básico de consumo gera equity premiums muito baixos (equity premium puzzle)
- R² muito baixo indica que o crescimento do consumo agregado explica pouca variação dos retornos
- Isto sugere que fatores além do consumo agregado são importantes para precificar ativos
- CAPM (Modelo de Um Fator):
- R² médio substancialmente maior, mostrando que o fator de mercado é relevante
- Consegue capturar parte significativa da variação cross-sectional dos retornos
- Porém, deixa espaço para fatores adicionais (como demonstrado pelo modelo FF3)
- Modelo Fama-French 3 Fatores:
- R² médio alto, demonstrando excelente poder explicativo
- Os fatores SMB (tamanho) e HML (valor) capturam dimensões importantes do risco não explicadas pelo mercado
- Confirma as anomalias empiricamente documentadas (small cap premium e value premium)
- Análise de Componentes Principais:
- Identifica fatores estatísticos que explicam a variação dos retornos sem impor estrutura econômica a priori
- Útil para comparar com fatores economicamente motivados (FF3)
- Fama-MacBeth:
- Testa se os prêmios de risco estimados são significativamente diferentes de zero
- Verifica se a relação risco-retorno é linear, conforme previsto pela teoria
Implicações:
- Os modelos multifatoriais (FF3) superam o CAPM e o CCAPM na explicação dos retornos
- Tamanho e valor são fatores de risco compensados no mercado
- O equity premium puzzle permanece um desafio teórico importante
Questão 2 - SARIMA
Principais Resultados:
- Estacionariedade:
- A série IPCA em nível não é estacionária (testes ADF e Phillips-Perron)
- A primeira diferença do logaritmo é estacionária, justificando d=1
- Forte componente sazonal identificado na FAC/FACP (picos em múltiplos de 12)
- Modelos Estimados:
- SARIMA(1,1,0)(1,0,0)₁₂: Componente sazonal autorregressivo
- SARIMA(1,1,0)(0,0,1)₁₂: Componente sazonal de médias móveis
- Ambos apresentam parâmetros significativos e resíduos sem autocorrelação (teste Ljung-Box)
- Raízes do polinômio característico fora do círculo unitário (estacionariedade/invertibilidade)
- Comparação AIC/BIC:
- Critérios de informação indicam qual modelo tem melhor balanço entre ajuste e parcimônia
- BIC tende a preferir modelos mais simples (penaliza mais a complexidade)
- Previsões Fora da Amostra:
- Rolling forecast implementado corretamente (reestimação passo a passo)
- RMSE permite comparar desempenho preditivo dos modelos
- O modelo com menor RMSE tem melhor capacidade preditiva
Implicações:
- Modelos SARIMA são apropriados para séries com sazonalidade (como IPCA)
- A diferenciação é necessária para obter estacionariedade
- A escolha entre SAR e SMA depende dos critérios de informação e desempenho preditivo
Referências:
- Cochrane, J. H. (2005). Asset Pricing. Princeton University Press.
- Fama, E. F., & French, K. R. (1993). Common risk factors in the returns on stocks and bonds. Journal of Financial Economics, 33(1), 3-56.
- Fama, E. F., & MacBeth, J. D. (1973). Risk, return, and equilibrium: Empirical tests. Journal of Political Economy, 81(3), 607-636.
- Sharpe, W. F. (1964). Capital asset prices: A theory of market equilibrium under conditions of risk. Journal of Finance, 19(3), 425-442.
- Mehra, R., & Prescott, E. C. (1985). The equity premium: A puzzle. Journal of Monetary Economics, 15(2), 145-161.
- Box, G. E. P., Jenkins, G. M., Reinsel, G. C., & Ljung, G. M. (2015). Time Series Analysis: Forecasting and Control. 5th ed. Wiley.
- Hyndman, R. J., & Athanasopoulos, G. (2021). Forecasting: Principles and Practice. 3rd ed. OTexts.
- Dickey, D. A., & Fuller, W. A. (1979). Distribution of the estimators for autoregressive time series with a unit root. Journal of the American Statistical Association, 74(366), 427-431.
- Phillips, P. C. B., & Perron, P. (1988). Testing for a unit root in time series regression. Biometrika, 75(2), 335-346.
Notas Metodológicas
Dados Utilizados:
- Portfólios Fama-French 25 (5×5): Disponíveis em Kenneth R. French Data Library
- Fatores Fama-French: Mkt-RF, SMB, HML do mesmo repositório
- Dados de Consumo: BEA NIPA Tables 2.1 (população), 2.3.4 (price index) e 2.3.5 (consumo nominal de não duráveis)
- Dados IPCA: IPEADATA (série número-índice mensal)
Pacotes R Utilizados:
tidyverse: Manipulação e visualização de dadosreadxl: Leitura de arquivos Excel (dados BEA)zoo/xts: Manipulação de séries temporaisforecast: Modelagem ARIMA/SARIMAtseries: Testes de raiz unitária (ADF, Phillips-Perron)knitr/kableExtra: Formatação de tabelas
Período de Análise:
- Questão 1 (Asset Pricing): 1947Q1 - 2025Q1 (dados trimestrais)
- Questão 2 (SARIMA): Janeiro/2000 - Janeiro/2025 (dados mensais)
Observações Importantes:
Os retornos dos portfólios Fama-French já são fornecidos em termos nominais (% ao mês). A conversão para trimestral usa agregação geométrica: \((1+r_1)(1+r_2)(1+r_3) - 1\).
Os retornos reais são calculados deflacionando pelos índices de preços do BEA, que já são agregados trimestralmente.
Todos os testes estatísticos utilizam nível de significância de 5% (α = 0.05).
As previsões out-of-sample na Questão 2 utilizam metodologia rolling window com reestimação a cada período.
Fontes de Dados
Questão 1 - Asset Pricing
Kenneth R. French Data Library: - 25 Portfolios Formed on Size and Investment (5×5): https://mba.tuck.dartmouth.edu/pages/faculty/ken.french/data_library.html - Fama/French 3 Factors: https://mba.tuck.dartmouth.edu/pages/faculty/ken.french/data_library.html
Bureau of Economic Analysis (BEA) - NIPA Tables: - Table 2.1: Personal Income and Its Disposition (população): https://apps.bea.gov/iTable/?reqid=19&step=2 - Table 2.3.4: Price Indexes for Personal Consumption Expenditures: https://apps.bea.gov/iTable/?reqid=19&step=2 - Table 2.3.5: Personal Consumption Expenditures by Major Type of Product: https://apps.bea.gov/iTable/?reqid=19&step=2
Questão 2 - Séries Temporais
IPEADATA: - IPCA - Índice Nacional de Preços ao Consumidor Amplo (número-índice mensal): http://www.ipeadata.gov.br/ - Código da série: PRECOS12_IPCAG12
Instruções de Reprodução
Pré-requisitos
Software: - R versão ≥ 4.0.0 - RStudio (recomendado) - Quarto CLI ≥ 1.3
Pacotes R necessários:
install.packages(c(
"tidyverse", "readxl", "zoo", "xts",
"forecast", "tseries", "knitr", "kableExtra",
"ggplot2", "gridExtra", "reshape2", "lubridate", "broom"
))Estrutura de Arquivos
Trabalho EPGE/
├── Lista2.qmd # Documento principal Quarto
├── Lista2.html # Saída renderizada
├── 25_Portfolios_5x5.csv # Dados dos 25 portfólios (French)
├── F-F_Research_Data_Factors.csv # Fatores Fama-French
├── BEA_Table_2.1.xlsx # População (BEA)
├── BEA_Table_2.3.4.xlsx # Price Index (BEA)
├── BEA_Table_2.3.5.xlsx # Consumo (BEA)
├── IPCA_IPEADATA.csv # Dados IPCA
└── README.md # Este arquivoPassos para Reprodução
- Download dos Dados:
- Baixe os arquivos CSV do Kenneth R. French Data Library
- Baixe as tabelas NIPA do BEA em formato Excel
- Baixe a série IPCA do IPEADATA em formato CSV
- Organização:
- Coloque todos os arquivos de dados no mesmo diretório do arquivo
Lista2.qmd - Verifique que os nomes dos arquivos correspondem aos especificados no código
- Coloque todos os arquivos de dados no mesmo diretório do arquivo
- Renderização:
- Abra
Lista2.qmdno RStudio - Execute:
quarto render Lista2.qmd - Ou use o botão “Render” no RStudio
- O arquivo HTML será gerado automaticamente
- Abra
- Tempo de Processamento:
- Primeira execução: ~5-10 minutos (cache vazio)
- Execuções subsequentes: ~2-3 minutos (com cache)
- O cache está habilitado para chunks computacionalmente intensivos
Resolução de Problemas
Erro de arquivo não encontrado: - Verifique os nomes dos arquivos de dados - Certifique-se de que estão no diretório correto - Use caminhos relativos ao arquivo .qmd
Erro de pacote não instalado: - Execute o comando de instalação dos pacotes acima - Reinicie a sessão R após instalação
Problemas com encoding: - Os arquivos CSV devem estar em UTF-8 - Use Sys.setlocale("LC_ALL", "en_US.UTF-8") se necessário
Informações de Sessão
Data de execução: 08/10/2025 16:43:54 Informações da Sessão R:R version 4.4.2 (2024-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64
Running under: Windows 11 x64 (build 26100)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.utf8
[2] LC_CTYPE=English_United States.utf8
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=English_United States.utf8
time zone: America/Sao_Paulo
tzcode source: internal
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] httr_1.4.7 jsonlite_1.8.9 urca_1.3-4 tseries_0.10-58
[5] forecast_8.23.0 reshape2_1.4.4 gridExtra_2.3 kableExtra_1.4.0
[9] knitr_1.49 broom_1.0.7 zoo_1.8-12 readxl_1.4.3
[13] lubridate_1.9.4 forcats_1.0.0 stringr_1.5.1 dplyr_1.1.4
[17] purrr_1.0.2 readr_2.1.5 tidyr_1.3.1 tibble_3.2.1
[21] ggplot2_3.5.2 tidyverse_2.0.0
loaded via a namespace (and not attached):
[1] gtable_0.3.6 xfun_0.52 htmlwidgets_1.6.4 lattice_0.22-6
[5] tzdb_0.4.0 quadprog_1.5-8 vctrs_0.6.5 tools_4.4.2
[9] generics_0.1.3 curl_6.1.0 parallel_4.4.2 xts_0.14.1
[13] pkgconfig_2.0.3 Matrix_1.7-1 lifecycle_1.0.4 compiler_4.4.2
[17] farver_2.1.2 munsell_0.5.1 codetools_0.2-20 htmltools_0.5.8.1
[21] yaml_2.3.10 pillar_1.10.1 nlme_3.1-166 fracdiff_1.5-3
[25] tidyselect_1.2.1 digest_0.6.37 stringi_1.8.4 labeling_0.4.3
[29] splines_4.4.2 fastmap_1.2.0 grid_4.4.2 colorspace_2.1-1
[33] cli_3.6.3 magrittr_2.0.3 withr_3.0.2 scales_1.3.0
[37] backports_1.5.0 timechange_0.3.0 TTR_0.24.4 rmarkdown_2.29
[41] quantmod_0.4.26 nnet_7.3-19 timeDate_4041.110 cellranger_1.1.0
[45] hms_1.1.3 evaluate_1.0.3 lmtest_0.9-40 viridisLite_0.4.2
[49] mgcv_1.9-1 rlang_1.1.4 Rcpp_1.0.13-1 glue_1.8.0
[53] xml2_1.3.6 svglite_2.1.3 rstudioapi_0.17.1 R6_2.6.1
[57] plyr_1.8.9 systemfonts_1.1.0Autores: José Victor
Disciplina: Econometria II - FGV EPGE
Data de Entrega: 10 de outubro de 2025