1 Load Libraries

# install.packages(c("tidyquant","timetk","lubridate","purrr",
#                    "quadprog","PerformanceAnalytics","ggplot2",
#                    "dplyr","tidyr","frenchdata","scales","knitr","kableExtra"))

library(tidyquant)
library(timetk)
library(lubridate)
library(purrr)
library(dplyr)
library(tidyr)
library(quadprog)
library(PerformanceAnalytics)
library(ggplot2)
library(frenchdata)
library(scales)
library(knitr)
library(kableExtra)

2 Q1 · Import ETF Daily Data

Download adjusted daily closing prices for 8 ETFs from Yahoo Finance (2010-01-01 to current date).

tickers <- c("SPY", "QQQ", "EEM", "IWM", "EFA", "TLT", "IYR", "GLD")

etf_raw <- tq_get(tickers,
                  from = "2010-01-01",
                  to   = Sys.Date(),
                  get  = "stock.prices")

# Wide tibble of adjusted prices
prices_tbl <- etf_raw %>%
  select(symbol, date, adjusted) %>%
  pivot_wider(names_from = symbol, values_from = adjusted) %>%
  arrange(date)

# Convert to xts for return calculations
prices_xts <- prices_tbl %>%
  tk_xts(date_var = date)

# First 6 rows
head(prices_xts) %>%
  as.data.frame() %>%
  round(2) %>%
  kable(caption = "Adjusted Prices — First 6 Trading Days") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE)

Adjusted Prices — First 6 Trading Days
	SPY	QQQ	EEM	IWM	EFA	TLT	IYR	GLD
2010-01-04	84.80	40.29	30.35	51.37	35.13	55.71	26.77	109.80
2010-01-05	85.02	40.29	30.57	51.19	35.16	56.07	26.83	109.70
2010-01-06	85.08	40.05	30.64	51.14	35.31	55.32	26.82	111.51
2010-01-07	85.44	40.07	30.46	51.52	35.17	55.41	27.06	110.82
2010-01-08	85.72	40.40	30.70	51.80	35.45	55.39	26.88	111.37
2010-01-11	85.84	40.24	30.64	51.59	35.74	55.08	27.01	112.85

# Last 6 rows
tail(prices_xts) %>%
  as.data.frame() %>%
  round(2) %>%
  kable(caption = "Adjusted Prices — Last 6 Trading Days") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE)

Adjusted Prices — Last 6 Trading Days
	SPY	QQQ	EEM	IWM	EFA	TLT	IYR	GLD
2026-06-02	759.57	746.16	70.80	291.66	105.02	85.65	99.99	411.95
2026-06-03	754.24	744.21	69.92	287.67	104.12	85.31	100.00	407.87
2026-06-04	757.09	740.61	69.10	292.01	104.95	85.50	101.79	411.27
2026-06-05	737.55	705.06	64.59	281.65	102.26	85.06	102.54	396.24
2026-06-08	739.22	716.07	65.75	284.11	102.88	84.62	101.08	397.27
2026-06-09	NA	NA	NA	NA	NA	NA	NA	NA

3 Q2 · Weekly and Monthly Simple Returns

# Weekly simple returns
weekly_returns <- prices_xts %>%
  apply.weekly(function(x) {
    r <- diff(x) / lag(x)
    tail(r, 1)
  }) %>%
  na.omit()

# Monthly simple returns
monthly_returns <- prices_xts %>%
  apply.monthly(function(x) {
    first_p <- as.numeric(x[1, ])
    last_p  <- as.numeric(x[nrow(x), ])
    matrix((last_p - first_p) / first_p,
           nrow = 1, dimnames = list(NULL, tickers))
  }) %>%
  na.omit()

cat("Weekly returns:  ", nrow(weekly_returns),  "observations\n")

## Weekly returns:   857 observations

cat("Monthly returns: ", nrow(monthly_returns), "observations\n")

## Monthly returns:  197 observations

tail(monthly_returns, 6) %>%
  as.data.frame() %>%
  mutate(across(everything(), ~ percent(., accuracy = 0.01))) %>%
  kable(caption = "Monthly Simple Returns — Last 6 Months") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE)

Monthly Simple Returns — Last 6 Months
	SPY	QQQ	EEM	IWM	EFA	TLT	IYR	GLD
2025-12-31	0.54%	-0.34%	2.23%	0.56%	3.22%	-1.43%	-1.09%	1.68%
2026-01-30	1.29%	1.43%	5.09%	4.37%	3.81%	0.11%	2.36%	11.72%
2026-02-27	-1.35%	-3.01%	5.57%	-0.29%	3.93%	4.93%	6.35%	13.26%
2026-03-31	-4.99%	-4.96%	-7.66%	-5.83%	-5.97%	-3.26%	-6.61%	-12.19%
2026-04-30	9.68%	14.28%	11.81%	11.38%	3.76%	-0.74%	8.17%	-3.23%
2026-05-29	4.97%	9.52%	6.97%	3.99%	2.64%	0.18%	-0.91%	-1.43%

4 Q3 · Monthly Returns as Tibble

monthly_tbl <- monthly_returns %>%
  tk_tbl(rename_index = "date") %>%
  mutate(date = as.yearmon(date))

head(monthly_tbl) %>%
  mutate(across(-date, ~ round(., 4))) %>%
  kable(caption = "Monthly Returns Tibble (first 6 rows)") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE)

Monthly Returns Tibble (first 6 rows)
date	SPY	QQQ	EEM	IWM	EFA	TLT	IYR	GLD
Jan 2010	-0.0524	-0.0782	-0.1037	-0.0605	-0.0749	0.0278	-0.0520	-0.0350
Feb 2010	0.0154	0.0347	-0.0089	0.0326	-0.0153	0.0057	0.0357	0.0100
Mar 2010	0.0500	0.0617	0.0631	0.0577	0.0556	-0.0201	0.0863	-0.0044
Apr 2010	0.0086	0.0224	-0.0271	0.0471	-0.0449	0.0358	0.0590	0.0463
May 2010	-0.0912	-0.0867	-0.0989	-0.0957	-0.1182	0.0525	-0.0852	0.0272
Jun 2010	-0.0355	-0.0510	0.0045	-0.0486	-0.0108	0.0506	-0.0278	0.0148

5 Q4 · Fama-French 3-Factor Data (Monthly)

ff3_raw <- download_french_data("Fama/French 3 Factors")

ff3_monthly <- ff3_raw$subsets$data[[1]] %>%
  rename(MktRF = `Mkt-RF`) %>%
  mutate(
    date  = as.yearmon(as.character(date), "%Y%m"),
    MktRF = MktRF / 100,
    SMB   = SMB   / 100,
    HML   = HML   / 100,
    RF    = RF    / 100
  ) %>%
  filter(date >= as.yearmon("2010-01"))

head(ff3_monthly) %>%
  mutate(across(-date, ~ round(., 5))) %>%
  kable(caption = "FF3 Factors — First 6 Months (decimal form)") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE)

FF3 Factors — First 6 Months (decimal form)
date	MktRF	SMB	HML	RF
Jan 2010	-0.0335	0.0043	0.0033	0e+00
Feb 2010	0.0339	0.0118	0.0318	0e+00
Mar 2010	0.0630	0.0146	0.0219	1e-04
Apr 2010	0.0199	0.0484	0.0296	1e-04
May 2010	-0.0790	0.0013	-0.0248	1e-04
Jun 2010	-0.0556	-0.0179	-0.0473	1e-04

6 Q5 · Merged Monthly Returns + FF3 Factors

merged_tbl <- monthly_tbl %>%
  inner_join(ff3_monthly, by = "date")

cat("Merged tibble:", nrow(merged_tbl), "rows,", ncol(merged_tbl), "columns\n")

## Merged tibble: 196 rows, 13 columns

cat("Date range:   ", as.character(min(merged_tbl$date)), "to",
                      as.character(max(merged_tbl$date)), "\n")

## Date range:    Jan 2010 to Apr 2026

head(merged_tbl) %>%
  mutate(across(-date, ~ round(., 4))) %>%
  kable(caption = "Merged Return + FF3 Tibble") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE, font_size = 11)

Merged Return + FF3 Tibble
date	SPY	QQQ	EEM	IWM	EFA	TLT	IYR	GLD	MktRF	SMB	HML	RF
Jan 2010	-0.0524	-0.0782	-0.1037	-0.0605	-0.0749	0.0278	-0.0520	-0.0350	-0.0335	0.0043	0.0033	0e+00
Feb 2010	0.0154	0.0347	-0.0089	0.0326	-0.0153	0.0057	0.0357	0.0100	0.0339	0.0118	0.0318	0e+00
Mar 2010	0.0500	0.0617	0.0631	0.0577	0.0556	-0.0201	0.0863	-0.0044	0.0630	0.0146	0.0219	1e-04
Apr 2010	0.0086	0.0224	-0.0271	0.0471	-0.0449	0.0358	0.0590	0.0463	0.0199	0.0484	0.0296	1e-04
May 2010	-0.0912	-0.0867	-0.0989	-0.0957	-0.1182	0.0525	-0.0852	0.0272	-0.0790	0.0013	-0.0248	1e-04
Jun 2010	-0.0355	-0.0510	0.0045	-0.0486	-0.0108	0.0506	-0.0278	0.0148	-0.0556	-0.0179	-0.0473	1e-04

7 Helper Functions

## ── Global Minimum Variance weights (long-only) ──────────────────────────────
gmv_weights <- function(cov_mat) {
  n    <- ncol(cov_mat)
  Dmat <- 2 * cov_mat
  dvec <- rep(0, n)
  # Equality: sum(w) = 1;  Inequality: w_i >= 0
  Amat <- cbind(rep(1, n), diag(n))
  bvec <- c(1, rep(0, n))
  sol  <- solve.QP(Dmat, dvec, Amat, bvec, meq = 1)
  w    <- sol$solution
  names(w) <- colnames(cov_mat)
  w
}

## ── CAPM residual covariance matrix ──────────────────────────────────────────
capm_cov <- function(ret_tbl, ff_tbl) {
  assets <- setdiff(names(ret_tbl), "date")
  df <- ret_tbl %>%
    inner_join(ff_tbl %>% select(date, MktRF, RF), by = "date") %>%
    mutate(across(all_of(assets), ~ . - RF, .names = "excess_{.col}"))

  betas  <- numeric(length(assets))
  resids <- matrix(NA, nrow(df), length(assets),
                   dimnames = list(NULL, assets))

  for (i in seq_along(assets)) {
    fit         <- lm(df[[paste0("excess_", assets[i])]] ~ df$MktRF)
    betas[i]    <- coef(fit)[2]
    resids[, i] <- residuals(fit)
  }

  sigma_m2   <- var(df$MktRF)
  beta_mat   <- matrix(betas, ncol = 1)
  cov_mat    <- (beta_mat %*% t(beta_mat)) * sigma_m2 + diag(apply(resids, 2, var))
  dimnames(cov_mat) <- list(assets, assets)
  cov_mat
}

## ── FF3 residual covariance matrix ───────────────────────────────────────────
ff3_cov <- function(ret_tbl, ff_tbl) {
  assets <- setdiff(names(ret_tbl), "date")
  df <- ret_tbl %>%
    inner_join(ff_tbl %>% select(date, MktRF, SMB, HML, RF), by = "date") %>%
    mutate(across(all_of(assets), ~ . - RF, .names = "excess_{.col}"))

  B      <- matrix(NA, 3, length(assets),
                   dimnames = list(c("MktRF","SMB","HML"), assets))
  resids <- matrix(NA, nrow(df), length(assets),
                   dimnames = list(NULL, assets))

  for (i in seq_along(assets)) {
    fit         <- lm(df[[paste0("excess_", assets[i])]] ~
                        df$MktRF + df$SMB + df$HML)
    B[, i]      <- coef(fit)[2:4]
    resids[, i] <- residuals(fit)
  }

  Sigma_F <- cov(df[, c("MktRF","SMB","HML")])
  cov_mat <- t(B) %*% Sigma_F %*% B + diag(apply(resids, 2, var))
  dimnames(cov_mat) <- list(assets, assets)
  cov_mat
}

8 Q6 · CAPM GMV — Covariance, Weights & Realized Return (2015/02)

Training window: 2010/02 – 2015/01 (60 months)
Allocation date: 2015/01 (weights computed)
Realized return: 2015/02

train_start  <- as.yearmon("2010-02")
train_end    <- as.yearmon("2015-01")
invest_date  <- as.yearmon("2015-02")

train_tbl <- merged_tbl %>%
  filter(date >= train_start & date <= train_end)

invest_ret <- merged_tbl %>%
  filter(date == invest_date) %>%
  select(all_of(tickers))

# CAPM covariance matrix
cov_capm_q6 <- capm_cov(
  train_tbl %>% select(date, all_of(tickers)),
  ff3_monthly
)

# GMV weights
w_capm_q6 <- gmv_weights(cov_capm_q6)

# Realized portfolio return
r_capm_q6 <- sum(w_capm_q6 * as.numeric(invest_ret))

# Weight table
data.frame(
  Asset  = names(w_capm_q6),
  Weight = percent(w_capm_q6, accuracy = 0.01)
) %>%
  kable(caption = "CAPM GMV Optimal Weights (as of 2015/01)") %>%
  kable_styling(bootstrap_options = c("striped","hover"),
                full_width = FALSE)

CAPM GMV Optimal Weights (as of 2015/01)
	Asset	Weight
SPY	SPY	26.04%
QQQ	QQQ	10.35%
EEM	EEM	0.52%
IWM	IWM	0.00%
EFA	EFA	3.48%
TLT	TLT	45.98%
IYR	IYR	5.78%
GLD	GLD	7.84%

# Weight bar chart
data.frame(Asset = names(w_capm_q6), Weight = w_capm_q6) %>%
  ggplot(aes(x = reorder(Asset, Weight), y = Weight, fill = Asset)) +
  geom_col(show.legend = FALSE) +
  geom_text(aes(label = percent(Weight, accuracy = 0.1)),
            hjust = -0.1, size = 3.5) +
  coord_flip() +
  scale_y_continuous(labels = percent_format(), limits = c(0, max(w_capm_q6) * 1.2)) +
  scale_fill_brewer(palette = "Set2") +
  labs(title = "Q6: CAPM GMV Weights (2015/01)",
       x = NULL, y = "Portfolio Weight") +
  theme_minimal(base_size = 12)

cat(sprintf("\n✅ Realized Portfolio Return (CAPM) in 2015/02: %.4f  (%.2f%%)\n",
            r_capm_q6, r_capm_q6 * 100))

## 
## ✅ Realized Portfolio Return (CAPM) in 2015/02: -0.0122  (-1.22%)

9 Q7 · FF3 GMV — Covariance, Weights & Realized Return (2015/02)

cov_ff3_q7 <- ff3_cov(
  train_tbl %>% select(date, all_of(tickers)),
  ff3_monthly
)

w_ff3_q7 <- gmv_weights(cov_ff3_q7)
r_ff3_q7 <- sum(w_ff3_q7 * as.numeric(invest_ret))

data.frame(
  Asset  = names(w_ff3_q7),
  Weight = percent(w_ff3_q7, accuracy = 0.01)
) %>%
  kable(caption = "FF3 GMV Optimal Weights (as of 2015/01)") %>%
  kable_styling(bootstrap_options = c("striped","hover"),
                full_width = FALSE)

FF3 GMV Optimal Weights (as of 2015/01)
	Asset	Weight
SPY	SPY	32.75%
QQQ	QQQ	2.82%
EEM	EEM	0.00%
IWM	IWM	2.92%
EFA	EFA	2.14%
TLT	TLT	46.89%
IYR	IYR	6.41%
GLD	GLD	6.08%

# Compare CAPM vs FF3 weights side by side
data.frame(
  Asset     = tickers,
  CAPM      = w_capm_q6,
  FF3       = w_ff3_q7
) %>%
  pivot_longer(-Asset, names_to = "Model", values_to = "Weight") %>%
  ggplot(aes(x = Asset, y = Weight, fill = Model)) +
  geom_col(position = "dodge") +
  scale_y_continuous(labels = percent_format()) +
  scale_fill_manual(values = c("CAPM" = "#2166ac", "FF3" = "#d6604d")) +
  labs(title = "Q6 vs Q7: CAPM and FF3 GMV Weights Comparison (2015/01)",
       x = NULL, y = "Portfolio Weight") +
  theme_minimal(base_size = 12) +
  theme(legend.position = "bottom")

cat(sprintf("\n✅ Realized Portfolio Return (FF3) in 2015/02:  %.4f  (%.2f%%)\n",
            r_ff3_q7, r_ff3_q7 * 100))

## 
## ✅ Realized Portfolio Return (FF3) in 2015/02:  -0.0132  (-1.32%)

# Side-by-side realized return summary
data.frame(
  Model             = c("CAPM GMV", "FF3 GMV"),
  Realized_Return   = percent(c(r_capm_q6, r_ff3_q7), accuracy = 0.01)
) %>%
  kable(caption = "Realized Portfolio Return in 2015/02") %>%
  kable_styling(bootstrap_options = c("striped","hover"),
                full_width = FALSE)

Realized Portfolio Return in 2015/02
Model	Realized_Return
CAPM GMV	-1.22%
FF3 GMV	-1.32%

10 Q8 · Rolling Backtest: CAPM vs FF3 GMV (2015/02 – 2026/05)

Strategy: At every month t, use 60-month history (t−60 to t−1) to estimate the GMV weights. Invest at t and record the realized return.

all_dates    <- sort(unique(merged_tbl$date))
roll_start   <- as.yearmon("2015-02")
roll_end     <- as.yearmon("2026-05")
invest_dates <- all_dates[all_dates >= roll_start & all_dates <= roll_end]

n_periods <- length(invest_dates)
port_returns_capm <- rep(NA_real_, n_periods)
port_returns_ff3  <- rep(NA_real_, n_periods)

for (i in seq_len(n_periods)) {
  t_invest      <- invest_dates[i]
  idx           <- which(all_dates == t_invest)
  t_end_train   <- all_dates[idx - 1]
  t_start_train <- all_dates[idx - 60]

  if (length(t_start_train) == 0 || is.na(t_start_train)) next

  win <- merged_tbl %>%
    filter(date >= t_start_train & date <= t_end_train)
  if (nrow(win) < 60) next

  ret_month <- merged_tbl %>%
    filter(date == t_invest) %>%
    select(all_of(tickers))
  if (nrow(ret_month) == 0) next

  tryCatch({
    cv_c <- capm_cov(win %>% select(date, all_of(tickers)), ff3_monthly)
    w_c  <- gmv_weights(cv_c)
    port_returns_capm[i] <- sum(w_c * as.numeric(ret_month))
  }, error = function(e) NULL)

  tryCatch({
    cv_f <- ff3_cov(win %>% select(date, all_of(tickers)), ff3_monthly)
    w_f  <- gmv_weights(cv_f)
    port_returns_ff3[i] <- sum(w_f * as.numeric(ret_month))
  }, error = function(e) NULL)
}

backtest_tbl <- tibble(
  date      = invest_dates,
  ret_capm  = port_returns_capm,
  ret_ff3   = port_returns_ff3
) %>%
  na.omit() %>%
  mutate(
    cum_capm = cumprod(1 + ret_capm) - 1,
    cum_ff3  = cumprod(1 + ret_ff3)  - 1,
    date_dt  = as.Date(date)
  )

cat("Backtest periods computed:", nrow(backtest_tbl), "\n")

## Backtest periods computed: 135

10.1 Cumulative Return Chart

backtest_tbl %>%
  select(date_dt, cum_capm, cum_ff3) %>%
  pivot_longer(-date_dt, names_to = "Model", values_to = "Cum_Return") %>%
  mutate(Model = recode(Model,
                        "cum_capm" = "CAPM GMV",
                        "cum_ff3"  = "FF3 GMV")) %>%
  ggplot(aes(x = date_dt, y = Cum_Return, color = Model)) +
  geom_line(size = 1.1) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "grey50") +
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  scale_color_manual(values = c("CAPM GMV" = "#2166ac",
                                "FF3 GMV"  = "#d6604d")) +
  labs(
    title    = "Cumulative Returns: CAPM GMV vs FF3 GMV",
    subtitle = "Rolling 60-month window · 8-Asset ETF Universe (2015/02 – 2026/05)",
    x        = NULL, y = "Cumulative Return", color = NULL,
    caption  = "Assets: SPY, QQQ, EEM, IWM, EFA, TLT, IYR, GLD"
  ) +
  theme_minimal(base_size = 13) +
  theme(plot.title    = element_text(face = "bold"),
        legend.position = "bottom")

10.2 Monthly Return Distribution

backtest_tbl %>%
  select(date_dt, ret_capm, ret_ff3) %>%
  pivot_longer(-date_dt, names_to = "Model", values_to = "Return") %>%
  mutate(Model = recode(Model,
                        "ret_capm" = "CAPM GMV",
                        "ret_ff3"  = "FF3 GMV")) %>%
  ggplot(aes(x = Return, fill = Model)) +
  geom_histogram(bins = 40, alpha = 0.7, position = "identity") +
  scale_x_continuous(labels = percent_format()) +
  scale_fill_manual(values = c("CAPM GMV" = "#2166ac",
                               "FF3 GMV"  = "#d6604d")) +
  labs(title = "Distribution of Monthly Portfolio Returns",
       x = "Monthly Return", y = "Count", fill = NULL) +
  theme_minimal(base_size = 12) +
  theme(legend.position = "bottom")

10.3 Performance Summary Table

perf_summary <- function(rets, label) {
  n   <- length(rets)
  ann_ret <- prod(1 + rets)^(12 / n) - 1
  ann_vol <- sd(rets) * sqrt(12)
  sr      <- ann_ret / ann_vol
  cum_ret <- cumprod(1 + rets) - 1
  drawdown <- cum_ret - cummax(cum_ret)
  max_dd  <- min(drawdown)
  tibble(
    Model              = label,
    `Ann. Return`      = percent(ann_ret, accuracy = 0.01),
    `Ann. Volatility`  = percent(ann_vol, accuracy = 0.01),
    `Sharpe Ratio`     = round(sr, 3),
    `Max Drawdown`     = percent(max_dd, accuracy = 0.01),
    `Total Months`     = n
  )
}

bind_rows(
  perf_summary(backtest_tbl$ret_capm, "CAPM GMV"),
  perf_summary(backtest_tbl$ret_ff3,  "FF3 GMV")
) %>%
  kable(caption = "Performance Summary: 2015/02 – 2026/05") %>%
  kable_styling(bootstrap_options = c("striped","hover","bordered"),
                full_width = FALSE) %>%
  row_spec(0, bold = TRUE)

Performance Summary: 2015/02 – 2026/05
Model	Ann. Return	Ann. Volatility	Sharpe Ratio	Max Drawdown	Total Months
CAPM GMV	5.50%	10.06%	0.546	-46.21%	135
FF3 GMV	5.30%	10.15%	0.522	-46.15%	135

End of Final Exam — Portfolio Analysis

Portfolio Analysis — Final Exam

Enkhmaa Gantulga

2026/06/10