# ============================================================
# Homework 04: GMVP & Tangency Portfolio for 4 Taiwan ETFs
# ETFs: 0050, 0056, 006205, 00646
# In-sample period: 2015/12/14 – 2018/12/28
# ============================================================
# install.packages(c("readr","quadprog")) # run once if needed
library(readr)
library(quadprog)
raw <- read_csv("myetf4.csv")
## Rows: 751 Columns: 5
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): Index
## dbl (4): tw0050, tw0056, tw006205, tw00646
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
colnames(raw) <- c("Date","ETF0050","ETF0056","ETF006205","ETF00646")
raw$Date <- as.Date(raw$Date, format = "%Y/%m/%d")
insample <- raw[raw$Date >= as.Date("2015-12-14") &
raw$Date <= as.Date("2018-12-28"), ]
cat("In-sample rows:", nrow(insample), "\n")
## In-sample rows: 751
# ============================================================
# Q1: GMVP using DAILY returns
# ============================================================
cat("\n============================================================\n")
##
## ============================================================
cat("Q1: GMVP based on DAILY returns\n")
## Q1: GMVP based on DAILY returns
cat("============================================================\n")
## ============================================================
prices_d <- as.matrix(insample[, c("ETF0050","ETF0056","ETF006205","ETF00646")])
# Compute daily log returns
ret_d <- log(prices_d[-1, ] / prices_d[-nrow(prices_d), ])
colnames(ret_d) <- c("0050","0056","006205","00646")
n <- ncol(ret_d)
mu_d <- colMeans(ret_d)
Sigma_d <- cov(ret_d)
cat("\nDaily mean returns:\n"); print(round(mu_d, 8))
##
## Daily mean returns:
## 0050 0056 006205 00646
## 0.00042381 0.00036185 -0.00027718 0.00022587
cat("\nDaily covariance matrix:\n"); print(round(Sigma_d, 8))
##
## Daily covariance matrix:
## 0050 0056 006205 00646
## 0050 7.894e-05 4.604e-05 0.00004514 3.696e-05
## 0056 4.604e-05 4.563e-05 0.00002695 2.375e-05
## 006205 4.514e-05 2.695e-05 0.00013103 2.939e-05
## 00646 3.696e-05 2.375e-05 0.00002939 5.914e-05
# --- GMVP via quadprog: min (1/2) w' Dmat w - dvec' w
# subject to: A' w >= bvec
# Enforce sum(w)=1 via two inequalities, long-only: w_i >= 0
Dmat <- 2 * Sigma_d
dvec <- rep(0, n)
# Each column of Amat is one constraint:
# col 1: +1'w >= 1 )
# col 2: -1'w >= -1 ) together: sum(w) = 1
# cols 3-6: I w >= 0 (long-only)
Amat <- cbind(rep(1,n), -rep(1,n), diag(n))
bvec <- c(1, -1, rep(0, n))
qp_d <- solve.QP(Dmat, dvec, Amat, bvec, meq = 0)
w_gmvp_d <- qp_d$solution
w_gmvp_d[w_gmvp_d < 1e-8] <- 0 # clean numerical noise
w_gmvp_d <- w_gmvp_d / sum(w_gmvp_d) # re-normalise
names(w_gmvp_d) <- colnames(ret_d)
cat("\n--- GMVP Weights (Daily, long-only) ---\n")
##
## --- GMVP Weights (Daily, long-only) ---
print(round(w_gmvp_d, 6))
## 0050 0056 006205 00646
## 0.000000 0.569975 0.083489 0.346536
cat("Sum of weights:", round(sum(w_gmvp_d), 6), "\n")
## Sum of weights: 1
ret_gmvp_d <- sum(w_gmvp_d * mu_d)
sd_gmvp_d <- sqrt(as.numeric(t(w_gmvp_d) %*% Sigma_d %*% w_gmvp_d))
cat(sprintf("\nGMVP Expected Daily Return : %.8f\n", ret_gmvp_d))
##
## GMVP Expected Daily Return : 0.00026138
cat(sprintf("GMVP Daily Std Dev : %.8f\n", sd_gmvp_d))
## GMVP Daily Std Dev : 0.00604032
# ============================================================
# Q2: GMVP using MONTHLY returns
# ============================================================
cat("\n============================================================\n")
##
## ============================================================
cat("Q2: GMVP based on MONTHLY returns\n")
## Q2: GMVP based on MONTHLY returns
cat("============================================================\n")
## ============================================================
# Extract last trading day of each calendar month
insample$YearMon <- format(insample$Date, "%Y-%m")
monthly_prices_df <- do.call(rbind,
lapply(split(insample, insample$YearMon), function(g) g[which.max(g$Date), ]))
monthly_prices_df <- monthly_prices_df[order(monthly_prices_df$Date), ]
prices_m <- as.matrix(monthly_prices_df[, c("ETF0050","ETF0056","ETF006205","ETF00646")])
ret_m <- log(prices_m[-1, ] / prices_m[-nrow(prices_m), ])
colnames(ret_m) <- c("0050","0056","006205","00646")
cat("Number of monthly return observations:", nrow(ret_m), "\n")
## Number of monthly return observations: 36
mu_m <- colMeans(ret_m)
Sigma_m <- cov(ret_m)
cat("\nMonthly mean returns:\n"); print(round(mu_m, 8))
##
## Monthly mean returns:
## 0050 0056 006205 00646
## 0.00820568 0.00662521 -0.00661152 0.00407547
cat("\nMonthly covariance matrix:\n"); print(round(Sigma_m, 8))
##
## Monthly covariance matrix:
## 0050 0056 006205 00646
## 0050 0.00119935 0.00087569 0.00086240 0.00040867
## 0056 0.00087569 0.00089960 0.00057302 0.00036511
## 006205 0.00086240 0.00057302 0.00260214 0.00070512
## 00646 0.00040867 0.00036511 0.00070512 0.00088550
Dmat_m <- 2 * Sigma_m
qp_m <- solve.QP(Dmat_m, dvec, Amat, bvec, meq = 0)
w_gmvp_m <- qp_m$solution
w_gmvp_m[w_gmvp_m < 1e-8] <- 0
w_gmvp_m <- w_gmvp_m / sum(w_gmvp_m)
names(w_gmvp_m) <- colnames(ret_m)
cat("\n--- GMVP Weights (Monthly, long-only) ---\n")
##
## --- GMVP Weights (Monthly, long-only) ---
print(round(w_gmvp_m, 6))
## 0050 0056 006205 00646
## 0.000000 0.493316 0.000000 0.506684
cat("Sum of weights:", round(sum(w_gmvp_m), 6), "\n")
## Sum of weights: 1
ret_gmvp_m <- sum(w_gmvp_m * mu_m)
sd_gmvp_m <- sqrt(as.numeric(t(w_gmvp_m) %*% Sigma_m %*% w_gmvp_m))
cat(sprintf("\nGMVP Expected Monthly Return : %.8f\n", ret_gmvp_m))
##
## GMVP Expected Monthly Return : 0.00533329
cat(sprintf("GMVP Monthly Std Dev : %.8f\n", sd_gmvp_m))
## GMVP Monthly Std Dev : 0.02507546
# ============================================================
# Q3: Tangency Portfolio (monthly, Rf = 0)
# ============================================================
cat("\n============================================================\n")
##
## ============================================================
cat("Q3: Tangency Portfolio (Monthly, Rf = 0)\n")
## Q3: Tangency Portfolio (Monthly, Rf = 0)
cat("============================================================\n")
## ============================================================
# Analytical formula (standard Markowitz; allows short-selling):
# w_tan proportional to Sigma^{-1} (mu - Rf * 1)
Rf <- 0
excess <- mu_m - Rf
Sigma_inv <- solve(Sigma_m)
z <- Sigma_inv %*% excess
w_tan <- as.numeric(z) / sum(z)
names(w_tan) <- colnames(ret_m)
cat("\n--- Tangency Portfolio Weights (Monthly, Rf=0) ---\n")
##
## --- Tangency Portfolio Weights (Monthly, Rf=0) ---
cat("(Unconstrained — short positions are allowed)\n")
## (Unconstrained — short positions are allowed)
print(round(w_tan, 6))
## 0050 0056 006205 00646
## 1.278969 -0.087371 -0.896758 0.705159
cat("Sum of weights:", round(sum(w_tan), 6), "\n")
## Sum of weights: 1
ret_tan <- sum(w_tan * mu_m)
sd_tan <- sqrt(as.numeric(t(w_tan) %*% Sigma_m %*% w_tan))
sharpe <- (ret_tan - Rf) / sd_tan
cat(sprintf("\nTangency Expected Monthly Return : %.8f\n", ret_tan))
##
## Tangency Expected Monthly Return : 0.01871875
cat(sprintf("Tangency Monthly Std Dev : %.8f\n", sd_tan))
## Tangency Monthly Std Dev : 0.04709414
cat(sprintf("Tangency Sharpe Ratio : %.8f\n", sharpe))
## Tangency Sharpe Ratio : 0.39747508
# ============================================================
# SUMMARY
# ============================================================
cat("\n============================================================\n")
##
## ============================================================
cat("SUMMARY OF RESULTS\n")
## SUMMARY OF RESULTS
cat("============================================================\n")
## ============================================================
cat("\n[Q1] GMVP (Daily Returns)\n")
##
## [Q1] GMVP (Daily Returns)
cat(" Weights:\n")
## Weights:
for(i in seq_along(w_gmvp_d))
cat(sprintf(" %-8s : %.6f\n", names(w_gmvp_d)[i], w_gmvp_d[i]))
## 0050 : 0.000000
## 0056 : 0.569975
## 006205 : 0.083489
## 00646 : 0.346536
cat(sprintf(" Portfolio Daily Return : %.8f\n", ret_gmvp_d))
## Portfolio Daily Return : 0.00026138
cat(sprintf(" Portfolio Daily Std Dev: %.8f\n", sd_gmvp_d))
## Portfolio Daily Std Dev: 0.00604032
cat("\n[Q2] GMVP (Monthly Returns)\n")
##
## [Q2] GMVP (Monthly Returns)
cat(" Weights:\n")
## Weights:
for(i in seq_along(w_gmvp_m))
cat(sprintf(" %-8s : %.6f\n", names(w_gmvp_m)[i], w_gmvp_m[i]))
## 0050 : 0.000000
## 0056 : 0.493316
## 006205 : 0.000000
## 00646 : 0.506684
cat(sprintf(" Portfolio Monthly Return : %.8f\n", ret_gmvp_m))
## Portfolio Monthly Return : 0.00533329
cat(sprintf(" Portfolio Monthly Std Dev: %.8f\n", sd_gmvp_m))
## Portfolio Monthly Std Dev: 0.02507546
cat("\n[Q3] Tangency Portfolio (Monthly, Rf=0)\n")
##
## [Q3] Tangency Portfolio (Monthly, Rf=0)
cat(" Weights (unconstrained):\n")
## Weights (unconstrained):
for(i in seq_along(w_tan))
cat(sprintf(" %-8s : %.6f\n", names(w_tan)[i], w_tan[i]))
## 0050 : 1.278969
## 0056 : -0.087371
## 006205 : -0.896758
## 00646 : 0.705159
cat(sprintf(" Portfolio Monthly Return : %.8f\n", ret_tan))
## Portfolio Monthly Return : 0.01871875
cat(sprintf(" Portfolio Monthly Std Dev: %.8f\n", sd_tan))
## Portfolio Monthly Std Dev: 0.04709414
cat(sprintf(" Sharpe Ratio : %.8f\n", sharpe))
## Sharpe Ratio : 0.39747508