1. Data Preparation

library(tidyverse)
library(lubridate)
library(quadprog)
library(knitr)
library(kableExtra)

# Load data
myetf4 <- read.csv("C:/Users/Admin/Downloads/myetf4.csv", stringsAsFactors = FALSE)

colnames(myetf4) <- c("Date", "tw0050", "tw0056", "tw006205", "tw00646")
myetf4$Date <- as.Date(myetf4$Date)

# Filter in-sample data
insample <- myetf4 %>%
  filter(Date >= as.Date("2015-12-14"),
         Date <= as.Date("2018-12-28")) %>%
  arrange(Date)

cat("In-sample rows:", nrow(insample), "\n")
## In-sample rows: 751
cat("Date range:", min(insample$Date), "to", max(insample$Date), "\n")
## Date range: 16783 to 17893
head(insample) %>%
  kable(caption = "First 6 rows of in-sample data") %>%
  kable_styling(bootstrap_options = c("striped", "hover"))
First 6 rows of in-sample data
Date tw0050 tw0056 tw006205 tw00646
2015-12-14 53.29 18.25 31.06 19.61
2015-12-15 53.33 18.38 31.59 19.63
2015-12-16 54.14 18.56 31.60 19.89
2015-12-17 54.77 18.81 32.23 20.05
2015-12-18 54.50 18.95 32.18 19.85
2015-12-21 54.41 19.02 33.00 19.64

2. Daily Returns

prices_daily <- insample[, c("tw0050", "tw0056", "tw006205", "tw00646")]

ret_daily <- as.data.frame(
  apply(prices_daily, 2, function(x) diff(x) / head(x, -1))
) * 100

summary(ret_daily) %>%
  kable(caption = "Summary of Daily Returns (%)") %>%
  kable_styling()
Summary of Daily Returns (%)
tw0050 tw0056 tw006205 tw00646
Min. :-7.03406 Min. :-5.43041 Min. :-6.83477 Min. :-4.10000
1st Qu.:-0.40630 1st Qu.:-0.28541 1st Qu.:-0.53456 1st Qu.:-0.30140
Median : 0.07183 Median : 0.08187 Median : 0.00000 Median : 0.04148
Mean : 0.04632 Mean : 0.03846 Mean :-0.02118 Mean : 0.02554
3rd Qu.: 0.50863 3rd Qu.: 0.41763 3rd Qu.: 0.48307 3rd Qu.: 0.41841
Max. : 3.73230 Max. : 2.48353 Max. : 6.00775 Max. : 5.63087

3. GMVP (Daily)

Model

Sigma_daily <- cov(ret_daily)

n <- ncol(ret_daily)

Dmat <- 2 * Sigma_daily
dvec <- rep(0, n)

Amat <- cbind(rep(1, n), diag(n))
bvec <- c(1, rep(0, n))

sol_daily <- solve.QP(Dmat, dvec, Amat, bvec, meq = 1)
w_daily <- sol_daily$solution
names(w_daily) <- colnames(ret_daily)

data.frame(ETF = names(w_daily),
           Weight = round(w_daily, 4)) %>%
  kable(caption = "GMVP Weights (Daily)")
GMVP Weights (Daily)
ETF Weight
tw0050 tw0050 0.0000
tw0056 tw0056 0.5718
tw006205 tw006205 0.0837
tw00646 tw00646 0.3445

Performance

mu_daily <- colMeans(ret_daily)

gmvp_ret_d <- sum(w_daily * mu_daily)
gmvp_sd_d  <- sqrt(t(w_daily) %*% Sigma_daily %*% w_daily)

data.frame(
  Metric = c("Expected Return (%)", "Std Dev (%)"),
  Value  = c(gmvp_ret_d, gmvp_sd_d)
) %>%
  kable(caption = "GMVP Performance (Daily)")
GMVP Performance (Daily)
Metric Value
Expected Return (%) 0.0290207
Std Dev (%) 0.6019007

4. Monthly Returns

insample_monthly <- insample %>%
  mutate(YM = format(Date, "%Y-%m")) %>%
  group_by(YM) %>%
  slice_tail(n = 1) %>%
  ungroup()

prices_monthly <- insample_monthly[, c("tw0050", "tw0056", "tw006205", "tw00646")]

ret_monthly <- as.data.frame(
  apply(prices_monthly, 2, function(x) diff(x) / head(x, -1))
) * 100

summary(ret_monthly) %>%
  kable(caption = "Summary of Monthly Returns (%)")
Summary of Monthly Returns (%)
tw0050 tw0056 tw006205 tw00646
Min. :-10.760 Min. :-7.6505 Min. :-17.30709 Min. :-8.7947
1st Qu.: -1.085 1st Qu.:-1.4821 1st Qu.: -2.97001 1st Qu.:-0.5192
Median : 1.195 Median : 0.5376 Median : -0.07805 Median : 1.1126
Mean : 0.882 Mean : 0.7087 Mean : -0.53555 Mean : 0.4511
3rd Qu.: 2.539 3rd Qu.: 2.3660 3rd Qu.: 2.48567 3rd Qu.: 2.3882
Max. : 6.047 Max. : 7.8421 Max. : 8.55715 Max. : 4.0964

5. GMVP (Monthly)

Sigma_monthly <- cov(ret_monthly)

Dmat <- 2 * Sigma_monthly
dvec <- rep(0, n)

Amat <- cbind(rep(1, n), diag(n))
bvec <- c(1, rep(0, n))

sol_monthly <- solve.QP(Dmat, dvec, Amat, bvec, meq = 1)
w_monthly <- sol_monthly$solution
names(w_monthly) <- colnames(ret_monthly)

data.frame(ETF = names(w_monthly),
           Weight = round(w_monthly, 4)) %>%
  kable(caption = "GMVP Weights (Monthly)")
GMVP Weights (Monthly)
ETF Weight
tw0050 tw0050 0.0032
tw0056 tw0056 0.4740
tw006205 tw006205 0.0012
tw00646 tw00646 0.5216

Performance

mu_monthly <- colMeans(ret_monthly)

gmvp_ret_m <- sum(w_monthly * mu_monthly)
gmvp_sd_m  <- sqrt(t(w_monthly) %*% Sigma_monthly %*% w_monthly)

data.frame(
  Metric = c("Expected Return (%)", "Std Dev (%)"),
  Value  = c(gmvp_ret_m, gmvp_sd_m)
) %>%
  kable(caption = "GMVP Performance (Monthly)")
GMVP Performance (Monthly)
Metric Value
Expected Return (%) 0.5733667
Std Dev (%) 2.4904413

6. Tangency Portfolio (Rf = 0)

Model

rf <- 0

Sigma_inv <- solve(Sigma_monthly)

z <- Sigma_inv %*% mu_monthly
w_tan <- z / sum(z)
names(w_tan) <- colnames(ret_monthly)

data.frame(ETF = names(w_tan),
           Weight = round(w_tan, 4)) %>%
  kable(caption = "Tangency Portfolio Weights")
Tangency Portfolio Weights
ETF Weight
tw0050 tw0050 1.3051
tw0056 tw0056 -0.1577
tw006205 tw006205 -0.8475
tw00646 tw00646 0.7002

Performance

tan_ret <- sum(w_tan * mu_monthly)
tan_sd  <- sqrt(t(w_tan) %*% Sigma_monthly %*% w_tan)
tan_sr  <- tan_ret / tan_sd

data.frame(
  Metric = c("Expected Return (%)", "Std Dev (%)", "Sharpe Ratio"),
  Value  = c(tan_ret, tan_sd, tan_sr)
) %>%
  kable(caption = "Tangency Portfolio Performance")
Tangency Portfolio Performance
Metric Value
Expected Return (%) 1.8090019
Std Dev (%) 4.4236381
Sharpe Ratio 0.4089399

7. Efficient Frontier

set.seed(123)

N <- 5000
sim_ret <- numeric(N)
sim_sd  <- numeric(N)

for (i in 1:N) {
  w <- runif(n)
  w <- w / sum(w)

  sim_ret[i] <- sum(w * mu_monthly)
  sim_sd[i]  <- sqrt(t(w) %*% Sigma_monthly %*% w)
}

plot(sim_sd, sim_ret,
     pch = 16, cex = 0.3, col = "lightblue",
     xlab = "Risk (Std Dev)",
     ylab = "Return",
     main = "Efficient Frontier")

points(gmvp_sd_m, gmvp_ret_m, col = "red", pch = 17, cex = 2)
points(tan_sd, tan_ret, col = "darkgreen", pch = 15, cex = 2)

legend("topleft",
       legend = c("Portfolios", "GMVP", "Tangency"),
       col = c("lightblue", "red", "darkgreen"),
       pch = c(16, 17, 15))